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ABSTRACT 


A  miniature  directional  sound  sensor  was  fabricated  using  micro-electro¬ 
mechanical  system  (MEMS)  technology  based  on  the  operational  principle  of 
Ormia  ochracea  fly’s  hearing  organ.  The  fly  uses  coupled  bars  hinged  at  the 
center  to  achieve  the  directional  sound  sensing  by  monitoring  the  difference  in 
vibration  amplitude  between  them.  The  MEMS  sensor  design  employed  in  this 
thesis  was  fabricated  using  the  PolyMUMPs  process.  The  sound  sensor  has  two 
primary  vibrational  modes  (rocking  and  bending)  which  were  simulated  by  finite 
element  analysis  and  tested  by  actuating  the  sensor  using  both  electrostatics  and 
sound.  The  experimentally  observed  vibrational  frequencies  were  found  to  be  in 
good  agreement  with  that  of  the  modeling  but  the  amplitudes  of  vibration  were 
found  to  be  relatively  small  compared  to  the  expected  values.  The  design  of 
optimized  sensors,  with  larger  amplitude  of  vibration  hence  high  sensitivity,  will 
be  discussed  in  this  presentation. 
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I.  INTRODUCTION 


“Conn,  Sonar,  possible  submerged  contact  bearing  030  or  290,”  bellows 
from  the  announcing  circuit  aboard  a  fast  attack  submarine.  “Helm,  right  ten 
degrees  rudder”,  the  Officer  of  the  Deck  orders  as  he  wonders  why  someone  had 
not  invented  a  better  way  of  pinpointing  which  side  of  the  towed  array  a  contact 
was  on.  Current  towed  array  technology  uses  an  array  of  omni-directional 
hydrophones.  Because  of  this,  the  towed  array  inherently  has  a  bearing 
ambiguity  that  can  only  be  resolved  using  the  Doppler  Effect  of  sound.  Turning 
the  submarine  either  toward  or  away  from  the  contact  would  cause  a  frequency 
shift  in  the  received  sound  and  thus  indicate  which  side  of  the  array  the  contact 
was  on.  If  frequency  went  up,  then  the  submarine  turned  towards,  and  if 
frequency  went  down,  it  turned  away.  There  are  two  major  problems  with  this 
technique  however.  First  it  assumes  that  the  contact  does  not  maneuver  while 
the  submarine  maneuvered,  which  is  not  always  the  case.  Secondly,  it  reduces 
the  submarine’s  stealth,  because  as  it  turns  through  the  water,  the  submarine 
itself  makes  noise.  If  a  directional  microphone  was  developed  and  perfected  it 
could  be  added  to  existing  towed  arrays,  and  could  eliminate  the  practice  of 
needing  to  turn  the  submarine  to  resolve  bearing  ambiguity,  thereby  improving 
submarine  stealth.  Eventually  new  arrays  could  be  built  completely  via  these 
microphones. 

Unmanned  Aerial  Vehicles  (UAV’s)  or  Unmanned  Underwater  Vehicles 
(UUVs)  currently  cannot  autonomously  turn  to  investigate  an  explosion  occurring 
on  its  right  or  left  side  for  further  analysis.  Existing  sensors  cannot  determine  the 
angle  the  sound  came  from  in  order  for  the  autonomous  system  to  decide  which 
way  to  turn  the  vehicle.  Unfortunately,  an  explosion  is  too  short  in  duration  to  use 
the  Doppler  Effect  to  determine  the  bearing  to  the  explosion;  so  current  searches 
must  involve  a  system  with  less  than  full  autonomy.  This  practice  could  be 
eliminated  if  an  array  of  directional  microphones,  small  enough  with  low  power 

consumption,  could  be  developed  and  perfected. 
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If  a  ground  soldier  is  being  fired  upon  by  a  sniper,  which  is  hidden  in  some 
aspect,  then  it  might  be  a  very  good  intelligence  point  if  the  bearing  of  the  gun 
fire  is  known.  If  a  directional  microphone  was  developed  and  perfected  it  could 
determine  the  bearing  and  thus  a  detailed  inspection  of  the  area  could  be 
conducted,  perhaps  saving  lives. 

There  are  a  multitude  of  places  in  the  military  where  it  might  be  important 
to  know  the  bearing  of  the  incident  sound.  Therefore  if  a  directional  microphone 
was  developed  and  perfected,  there  are  endless  applications  within  the  military. 
Basically,  anywhere  there  is  a  need  to  know  what  direction  the  sound  came  from, 
on  a  low  cost,  small  scale,  and  low  power  consuming  aspect  this  perfected 
directional  microphone  is  the  military’s  solution. 

A.  BACKGROUND 

As  is  often  the  case,  biology  has  solved  the  sound  bearing  problem  long 
ago.  It  turns  out  that  “Mother  Nature”  has  developed  a  biological  system  that  can 
not  only  hear  the  sound,  but  can  determine  the  angle  from  which  the  sound 
came. 


1.  Biological  Motivation 

There  exists  a  fly,  the  Ormia  ochracea,  which  has  ears  that  are  extremely 
sensitive  to  the  direction  of  sound  and  are  small  enough  that  they  can  easily  be 
mimicked  within  a  biomimetic  silicon  chip. 

The  genus  Ormia  is  a  parasitoid  and  procreates  by  the  female  fly  laying  its 
larvae  on  live  male  crickets  [Miles  et  at.,  1995],  Most  common  hosts  are  genus 
Gryiius  or  Scapteriscus;  for  example  in  the  case  of  Ormia  ochracea,  it  is  male 
field  crickets  [e.g.  Gryiius  firmus]  that  are  the  hosts  [Christe,  2003],  while  with 
Ormia  depleta,  it  is  the  male  mole  crickets  [e.g.  Scapteriscus  borellii ]  that  are  the 
hosts  [Frank  et  at.,  2006],  The  larvae  then  eat  the  cricket  for  nourishment  in 
preparation  for  their  metamorphosis  [Robert  et  ai,  1994],  However,  crickets  are 
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mostly  nocturnal  and  black  so  the  fly  cannot  rely  on  its  eyesight  to  find  them. 
Instead  it  has  developed  ears  that  not  only  hear  the  male  crickets  mating  call,  but 
also  tell  the  fly  the  angle  the  sound  came  from,  so  the  fly  can  vector  in  on  the 
location  of  the  cricket.  For  ease  of  understanding  the  mechanics  of  the  biology, 
only  the  species  Ormia  ochracea’s  ears  will  be  explained,  since  it  will  otherwise 
become  too  confusing  to  discuss  the  minute  differences  between  all  the  Ormia 
species’  ears. 


a.  Fly’s  Auditory  Mechanics 

The  fly  itself  is  relatively  small,  only  measuring  about  a  centimeter 
in  length  [Christe,  2003],  Because  of  its  small  size,  the  ear  structure  is  not  very 
big  either;  the  ear  drums  are  separated  by  a  mere  500  pm  [Christe,  2003],  The 
ears  of  the  fly  are  not  located  on  the  head  as  in  humans,  but  instead  they  are 
located  on  the  thorax  just  behind  the  head  [Miles  et  al. ,  1995],  Their  location 
along  with  their  small  size  means  that  the  fly  cannot  rely  on  the  difference  in 
arrival  times  and  amplitudes  of  sound  at  each  ear,  as  most  animals  do,  including 
humans.  Instead  the  fly  developed  a  unique  way  of  amplifying  the  arrival  time 
difference  and  amplitude  difference.  The  fly  achieves  this  difference  by 
mechanically  coupling  both  of  its  ear  drums  together  via  a  bridging  membrane 
[Miles  et  al.,  1995], 

For  those  that  might  enjoy  the  biological  terms  and  all  those 
scientific  names,  specific  details  of  the  internals  of  the  ear  structure  can  be  read 
about  in  Robert  et  al.  (1994).  However,  for  physics  majors,  a  more  basic 
description  of  the  fly’s  ears  is  more  appropriate.  Looking  at  the  fly’s  auditory 
organ,  as  in  Figure  1  below,  the  ear  structure  can  be  broken  down  into  two 
mechanical  bars  connected  by  a  hinge  membrane  [Miles  et  al.,  1995],  Parts  1 
and  2  from  the  picture  are  the  bars,  part  3  is  the  hinge.  PTM  is  prosternal 
tympanal  membrane  and  for  all  intents  and  purposes  is  the  ear  drum  of  the  fly 
and  is  the  mechanism  that  drives  the  solid  bars  via  the  sound  input  [Miles  et  al., 
1995], 
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Figure  1 .  Fly’s  auditory  system  and  mechanical  components  (After:  Miles  et 
al.,  1995). 

The  key  physics  behind  how  the  ear  structure  works  so  well  at 
determining  the  angle  of  incidence  is  the  solid  bar  and  hinge  interaction.  This 
interaction  can  be  modeled  via  a  simplified  yet  accurate  mechanical  equivalent 
system. 


b.  Simplified  Mechanical  System 

A  quick  glance  at  this  mechanical  system  makes  it  seem  difficult 
and  overwhelming,  but  the  mechanical  system  can  be  simplified  into  a  simple 
schematic  of  springs,  bars,  damping  devices,  etc.  as  shown  below  in  Figure  2 
[Miles  et  al.,  1995], 
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Figure  2.  Simplified  mechanical  equivalent  of  the  fly’s  ear  system  (From:  Miles 
etai,  1995). 

This  is  essentially  a  two  degree  of  freedom  spring-mass-damper 
system  with  a  torsional  element.  To  fully  understand  the  fly’s  mechanism,  a  more 
detailed  look  into  two  degree  of  freedom  spring-mass-damper  systems  is  now 
appropriate. 

B.  PHYSICS  OF  WORKING  SYSTEM 

1.  Simple  Two  Degrees  of  Freedom  Spring-Mass-Damper  System 

Assume  for  ease  of  understanding  that  the  mechanical  system  of  interest 
is  modeled  as  in  Figure  3,  a  simple  two  degree  of  freedom  spring-mass-damper 
system. 
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Figure  3.  Two  degrees  of  freedom  spring-mass-damper  system  (After:  Rao, 
2003). 

Then  the  free-body  diagrams  look  like  those  in  Figure  4. 


Spring  ki  under  tension  for  +  X!  Spring  k3  under  tension  for  +  (x^)  Spring  k2  under  compression  for  +  x2 

Figure  4.  Free-body  diagram  (After:  Rao,  2003). 

From  these  diagrams  the  equations  of  motion  can  be  determined  applying 
Newton’s  second  law  as  discussed  by  S.  Rao  [Rao,  2003], 

mx  =  Fnel  (1.1) 
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Think  about  all  the  forces  on  each  body  in  Figure  4  and  place  them  on  the 
right  side  of  the  equation.  The  force  due  to  a  spring  is  opposite  in  direction  and 
proportional  to  the  displacement;  while  the  damping  force  is  opposite  in  direction 
and  proportional  to  the  speed.  Therefore  use  the  following  equations  to  define 
the  forces. 


F 


spring 


-he;  F 


damping 


-C  x 


(1.2) 


The  k  is  known  as  the  spring  constant  and  the  C  is  known  as  the  damping 
constant.  Thus,  the  equations  of  motion  for  the  system  in  Figure  4  become: 

mxx j  =-k1xl  -C,i,  +  k3(x2-x1)  +  C3(x2-xl)  +  Fl  (1 .3) 

^2^2  ~  —^3 (x2  —  x^  —  k3 ( x,  —  x^  —  k2x2  —  C2x2  +  F2  (14) 

It  is  important  to  note  that  if  the  system  is  moved  in  the  positive  X  direction  then 
spring  k2  will  be  under  compression  and  therefore  a  negative  force  too.  Hence 
the  two  damping  forces  and  the  two  spring  forces  are  negative  for  equation  1 .4, 
vice  the  two  positive  and  two  negative  in  equation  1.3,  where  the  springs  were 
both  under  tension. 


Grouping  like  terms  and  rearranging,  these  can  be  rewritten  as  [Rao, 


2003]: 


mxxx  +  (Cx  +  C3)xx  -  C2x2  +  (k}  +  ki)xl  - k3x2  =  Fx 
m2x2  - C3x j  +  (C3  +  C2)x2  - kyX}  +  (k,  +  k2)x2  =  F2 


(1.5) 

(1.6) 


These  are  two  coupled  second  order  ordinary  differential  equations  and 
can  be  expressed  in  matrix  form  as  [Rao,  2003]: 

[, m]x(t )  +  [C]i(?)  +  [A']T(?)  =  F(t )  (1 .7) 

Where  the  mass,  damping,  and  spring  constant  (stiffness)  matrices  are: 

[' m ]  ■■ 


f  m  0  ^ 

(' c.  +  c, 

-c3  ) 

k  k.+k,  -kx  \ 

I 

;  [C]  = 

1  J 

;  [k]  = 

1  3  3 

l  0  mi) 

l  ~C3 

c3  +  cj 

While  the  displacement  and  force  matrices  are: 
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^i(0^ 

V*2(0j 

;  F(t)  = 

U(0, 

x{t)  = 

The  displacements  of  the  two  masses  are  in  general  expressed  as  a  linear 
combination  of  the  Eigen  modes  [Rao,  2003], 


The  Eigen  modes  are  calculated  by  solving  the  differential  equation  1.7. 
To  simplify  the  math  and  to  demonstrate  the  process,  assume  there  is  no 
damping  in  equation  1.7.  Doing  so,  equation  1.7  becomes  equation  1.8. 

[m]x(t)  +  [k]x(t)  -  F(t)  (1.8) 


Solving  the  homogenous  solution  the  equation  simplifies  to  equation  1 .9. 

[m]x(t)  +  [k]x(t)  =  0  (1.9) 


Assuming  the  solution  takes  the  form  of  aeiM  as  is  done  in  Taylor’s  classical 
mechanics  book,  then  equation  1.9  becomes  equation  1.10  [Taylor,  2005], 

-o)2[m]aeiM  =-[k]aeirot  (1.10) 

Cancel  like  terms  from  both  sides  and  rearrange  to  take  the  form  of  equation 
1.11  [Taylor,  2005], 

([k]-o)2[m])a  =  0  (1.11) 


At  this  point  it  is  apparent  that  “a”  could  be  zero,  but  that  is  a  trivial  solution.  So 
the  next  logical  solution  is  that  the  determinant  of  ([k]-co2[m])  must  be  zero 


[Taylor,  2005],  The  determinant  of  ([£]- co2[m\)  is  ( k-mor^[k  +  k3-mco 2) 

assuming  ki=k2=k  [Taylor,  2005],  Therefore  the  first  two  natural  frequencies  or 
Eigen  modes  are: 


<x>, 


co2 


k  +  2k3 


m 


(1.12) 


In  the  following  section,  the  response  of  the  Fly’s  ears  will  be  analyzed 
using  the  above  matrix  formalism. 
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2.  Fly’s  Mechanical  Model 


Just  as  was  done  for  the  simple  two  degree  of  freedom  example,  the  first 
step  is  to  take  the  mechanical  system  of  the  fly  in  Figure  5  and  produce  an 
equivalent  free-body  diagram.  Since  the  right  and  left  ear  structures  are  the 
same,  Ki=K2  and  Ci=C2;  it  is  just  easier  to  call  them  Ks  and  Cs  (s  for  structure).  In 
addition  k3  and  c3  where  renamed  to  Kt,  Ct  due  to  the  coupling  membrane  (t  for 
tympanal  bridge).  Notice  that  the  spring  Kt  expands  as  and  springs  Ks  get 
compressed  as  0i  and  02  get  bigger.  This  is  an  important  difference  from  the 
simple  two  degree  of  freedom  example  in  Figure  3. 


Kt,  Ct 


Figure  5.  Simplified  mechanical  equivalent  of  the  fly’s  hearing  System  (After: 

Miles  et  al.,  1995). 

The  amount  of  displacement  of  the  fly’s  ears  is  very  small  (nanometer) 
compared  to  the  length  (L)  of  the  bar  (millimeter);  therefore  the  small  angle 
approximation  applies.  That  is  instead  of  X=L*sin(0),  X  can  be  approximated  to 
be  X  =  L*0  and  the  analysis  can  be  done  in  terms  of  X  vice  0.  Keeping  the 
generalized  coordinate  as  X,  vice  shifting  to  0,  means  that  the  fly  system  can  be 
more  easily  compared  to  the  simple  two  degrees  of  freedom  example  discussed 
earlier. 
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Figure  6  is  the  free-body  diagram  produced  from  the  mechanical  system  in 
Figure  5.  Springs  Ks  where  compressed,  due  to  F i  and  F2,  and  thus  have  a 
restoring  force  upward.  Spring  Kt  causes  a  restoring  force  upward  on  each  bar  as 
well.  The  damping  forces  too  are  upward. 


KsX, 


Cs*! 


K,(Xi+X2) 


csx2 


KsX? 


Figure  6.  Fly  free-body  diagram. 

The  next  step  is  to  formulate  the  equations  of  motion  for  the  system  using 


Newton’s  second  law  again  as  follows: 

=  ~Ksx i - CA  - K, (xi  +xi)~Ct{x2+xl)  +  Fl  (1.13) 

m2x2  =  -Ksx2  -  Csx2  -KXx2+xx)-CXx2+xx)  +  F2  (1.1 4) 

Rearranging  and  grouping  like  terms  as  done  before  in  the  simple  two  degree  of 
freedom  example,  the  equations  become: 

mlxl  +  (Cs  +  Ct  )x,  +Ctx2+  (Ks  +Kt)xl+Ktx2=Fl  (1-15) 

w2x 2  +  (Cs  +  Ct )  +  CrXj  +  ( Ks  +  Kt  )x2  +  Ktxx  =  F2  (1 .16) 


Since  the  bars  are  identical  their  masses  are  equal,  m-i=m2,  and  thus  the  coupled 
equations  of  motion  in  matrix  form  becomes  equation  1.12  [Miles  et  al.,  1995], 
Equation  1 . 1 7  is  the  primary  equation  for  the  motion  of  the  fly’s  ear  structure. 


m 

O' 

"xi" 

+ 

~c,+ct 

c, 

"Xj" 

+ 

~KS+Kt 

Kt 

"xi" 

>r 

0 

m 

x2 

Ct 

Cs+Ct_ 

x2 

K, 

Ks  +  Kt_ 

x2 

i 

to  J 

1 _ 

10 


The  torsional  form  of  this  equation  is  not  as  straight  forward  but  included 
here  as  background  for  issues  relating  to  inertia  to  be  discussed  later.  The 
transformation  from  equation  1.17  to  1.19  used  the  following  equations: 

X  =  L6 

F  =  T/L  (1.18) 

Ism  =  / 


By  placing  the  two  equations  (1.17  and  1.19)  close  together  it  is  easier  to  see  the 
translation  from  mass  to  inertia  and  force  to  torque. 


Inertia  Damping  Stiffness 

I  Olpil  2r  C,+Ct  Ct  W  J  Ks+Kt  Kt 

0  l\[a J  [  Ct  Cs+Ct\[02\  [  Kt  Ks+Kt 


Torque 


3.  Eigen  Frequency  Modes 


As  was  discussed  for  the  simple  two  degree  of  freedom  system,  the  fly  too 
has  a  first  and  second  mode  of  its  natural  frequencies.  These  are  called  the 
rocking  mode  and  the  bending  mode.  The  rocking  mode  is  called  that  because 
the  motion  of  the  coupled  bar  mechanism  looks  like  a  seesaw  as  it  rocks  back 
and  forth.  The  bending  mode  is  so  called  because  it  looks  like  a  bird  bending  its 
wings  or  flapping  its  wings  for  flight.  Look  at  Figure  7  for  a  visual  of  the  bars  at  an 
instance  in  time  for  each  mode. 
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MODE 

Figure  7.  Ear  vibration  Eigen  modes  (After:  Robert  et  at.,  1996). 


a.  Rocking  Mode  and  Bending  Modes 


The  Eigen  modes  of  the  coupled  system  can  be  obtained  following 
the  process  as  described  in  the  two  degrees  of  freedom  section  but  not 
neglecting  damping.  The  exact  solution  can  be  obtained  by  solving  equation  1.17 
as  described  in  the  Miles  et  al  paper.  The  first  two  modes  are  as  follows  [Miles  et 
at.,  1995]: 


^Rocking 


(1.20) 


co 


Bending 


K+2K. 


M 


(1.21) 


The  rocking  motion  depends  only  on  the  two  springs  at  the  edges  in 
Figure  5,  while  the  bending  motion  involves  springs  at  the  edges  as  well  as  at  the 
center. 


4.  Displacement  due  to  Sound 

Which  mode  is  activated  during  oscillations  of  the  structure  is  purely 
dependent  on  the  initial  conditions.  For  example,  if  the  sound  is  incident  from 
directly  above  the  structure  where  the  angle  of  incidence  is  zero,  only  the 
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bending  mode  can  be  actuated.  However,  the  most  probable  situation  is  that 
initial  conditions  are  such  that  linear  combinations  of  the  two  modes  occur 
simultaneously. 

If  the  sound  is  incident  from  the  left  at  an  angle  as  shown  in  Figure  7,  the 
left  ear  structure  is  undergoing  a  linear  combination  of  the  rocking  and  bending 
modes.  The  displacements  due  to  these  two  modes  add  together  to  create  a 
large  amplitude  on  the  left  side.  While  for  the  right  ear  structure,  the  rocking 
mode  pushes  the  right  bar  up  and  the  bending  mode  moves  it  down.  These  two 
subtract  from  each  other  creating  small  amplitude  on  the  right  ear.  The  net  result 
is  a  difference  in  amplitudes  depending  on  the  angle  of  incidence  of  the  sound. 


Expressing  this  in  equation  form  the  linear  combination  becomes 
equations  1.22  and  1.23  for  the  left  (x-i)  and  right  (x2)  respectively  [Miles  et  al., 
1995],  The  subscript  b  is  for  bending  mode,  r  is  for  rocking  mode,  and  co  is  the 
excitation  frequency. 


The  amplitudes  of  the  displacement  for  each  mode  are  expressed  in 
equations  1.24  and  1.25  [Miles  et  al.,  1995];  m  and  s  are  the  mass  and  area  of 
the  ear  structure,  t  is  the  time  delay  between  the  sound  force  on  the  left  and  right 
bar,  and  finally  P  is  the  amplitude  of  sound  pressure  which  hardly  varies  across 
the  fly’s  ear.  The  specific  derivation  was  given  by  Miles  et  al,  1995.  Modifying 
their  equations  the  phase  constants  for  rocking  and  bending  modes  can  be 
expressed  in  terms  of  a  damping  coefficient  y. 

xl  (f)  =  Ab  sin  (cot +  <f>h)  +  Ar  cos(cot  +  <f>r) 
x2{t)  =  Ab  sin(W  +  ^)  -  4.c°s(^  +  ^) 


(1.22) 

(1.23) 


A.  =  — 


Ps 

m 


(  \ 

sin  (&>r/2) 

^j(co;  -co2f  +(r,.co)2 


(1.24) 
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A 


b 


\ 


(1.25) 


(tr) 


f 


=  arc  tan 


y,.co 


2  2 


(K) 


f 


=  arc  tan 


Zb® 


(col -CO2) 


(1.26) 

(1.27) 


Typical  values  of  the  parameters  for  this  mechanical  model  of  the  fly  come 
from  the  Miles  et  al. ,  1995,  paper.  Using  these  values,  the  amplitudes  of 
oscillation  for  the  right  and  left  ears  based  on  equations  1.21  and  1.22  were 
calculated  using  MATLAB.  Figure  8  shows  amplitudes  and  phases  of  the  left  and 
right  ears  as  a  function  of  sound  frequency.  In  this  analysis  sound  is  assumed  to 
be  incident  at  45  degrees  with  a  pressure  of  one  Pascal. 

For  example,  a  5  kHz  sound  wave,  which  arrives  at  both  ears  with  the 
same  magnitude,  is  detected  with  approximately  a  10  dB  magnitude  difference 
by  the  fly’s  ears  [Karunasiri  et  al.,  2005], 
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Figure  8.  Frequency  response  of  the  ears  based  on  the  model  (From: 
Karunasiri  et  at.,  2005). 

Figure  9  shows  the  time  dependence  of  oscillations  with  relatively  large 
amplitude  difference  created  when  the  incident  sound  comes  from  a  45  degree 
angle.  The  phase  shift  is  clearly  visible  from  the  two  traces  in  Figure  9. 


15 


Figure  9.  Transient  response  based  on  the  model  for  45  degree  inclined  5kHz 
sound  wave  (After:  Karunasiri  et  a!.,  2005). 


Analytical  expressions  of  the  two  amplitudes  oscillation  in  Figure  9  can  be 
obtained  using  the  two  Eigen  modes.  The  displacements  of  the  two  sides  can  be 
written  as: 


x,(r)  =  A1  sin^cot  +  fa) 


(1.28) 


x2(f)  =  A2  sin(<xtf  +  02) 


(1.29) 


where  A-i  and  A2  are  the  amplitudes  of  the  left  and  right  ears  and  </>  i  and  <j)2  are 
the  corresponding  phases.  Using  equations  1 .22  to  1 .29,  it  can  be  shown  that  the 
two  amplitudes  and  phase  difference  are  given  by: 
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Al  +  A;+2AhArsm{<t>h-<t>r) 


(1.30) 


4  =  4  +  4  -  2 44  sin  ( <Ph  -  </>r ) 


cos  (fc-fa) 


4-4 

44 


(1.31) 

(1.32) 


It  can  be  easily  seen  that  the  amplitude  ratio  and  phase  difference  are 
independent  of  the  incident  sound  pressure. 

5.  Determination  of  the  Direction  of  Sound 

It  has  been  shown  that  the  fly’s  mechanical  system  responds  with  a 
difference  in  amplitude  for  the  side  that  the  sound  is  incident  from,  it  is  also 
important  to  discuss  that  there  is  an  amplification  in  the  sensed  time  difference  of 
arrival  of  the  sound  waves  too. 

a.  Sensed  Time  Difference 

If  the  fly  did  not  amplify  this  signal,  the  time  arrival  difference  would 
only  be  approximately  2.5  micro  seconds  [Karunasiri  et  a!.,  2005],  This  is  too 
small  for  the  fly  to  notice  an  appreciable  difference.  However,  according  to  Figure 
9,  there  is  a  large  phase  difference  between  the  left  and  right  ears.  This  phase 
difference  corresponds  to  a  sense  time  difference  of  about  50  microseconds 
[Karunasiri  et  at.,  2005],  That  corresponds  to  an  amplification  of  20  times.  Figure 
10  shows  the  amplification  of  the  sensed  time  difference  for  different  frequencies. 

What  is  most  noticeable  is  the  linear  region  from  -40  to  +40 
degrees  for  5  kHz  tonal.  It  is  no  coincidence  that  biology  has  evolved  the 
parasitic  fly’s  ears  to  the  tonal  frequency  for  the  cricket  host.  The  fly  shows 
excellent  directionality  for  this  region  of  incident  angles  and  can  thus  vector  in 
easily  on  its  prey. 
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So  if  the  frequency  of  the  incoming  signal  is  known,  and  the  time 
difference  is  sensed,  then  it  is  possible  to  derive  the  angle  of  incident  sound 
using  equation  1.33,  with  relative  ease. 


(^2  h ) 


—  cos 

CO 


42-42 

AxA2 


(1.33) 


Figure  10.  Sensed  time  difference  for  different  frequencies  and  angles. 

C.  MATLAB  PROGRAM  TO  VISUALIZE  SYSTEM  WORKING 

In  order  to  visualize  the  fly’s  bars  in  action,  a  MATLAB  simulation  program 
was  created  using  the  equations  from  the  Miles  et  al.,  1995,  paper  and  the 
parameter  values  from  the  Roberts  et  al.,  1994,  paper.  Code  for  the  MATLAB 
program  is  contained  in  Appendix  A,  so  that  the  program  can  be  run  and  seen 

working  with  the  ear  bars  actually  moving,  if  so  desired. 
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1.  Sound  from  Zero  Degrees 


Figures  10,  11,  and  12  show  the  bars  at  chronological  instances  in  time  for 
sound  incident  from  zero  degrees.  They  are  actual  screen  shots  while  the 
program  was  running  in  a  PowerPoint  presentation.  Notice  that  there  is  no 
amplitude  difference  between  the  left  (blue)  and  right  bars  (red). 


Figure  1 1 .  Zero  degrees  incidence  with  bars  up. 


Figure  12.  Still  zero  degrees  and  the  bars  moved  down. 
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Figure  13.  Zero  degrees  and  bars  moved  completely  down. 

The  motion  depicted  in  the  transition  from  Figure  1 1  to  Figure  13  is  the 
bending  mode  and  the  animation  looks  like  a  bird  flapping  its  wings.  Again  notice 
that  as  the  wings  flapped,  the  amplitude  on  each  side  was  the  same. 

2.  Sound  from  30  Degrees 

Figures  14,  15,  and  16  show  the  same  chronological  transition  in  time  but 
for  sound  incident  from  30  degrees.  Here  the  motion  will  not  be  pure  bending 
mode,  but  a  linear  combination  of  the  rocking  mode  and  the  bending  mode.  It 
looks  like  the  wings  are  flapping,  but  the  right  wing  flaps  with  a  higher  amplitude 
and  phase. 

Note  that  in  this  MATLAB  program  the  sound  is  coming  from  the  right, 
unlike  diagram  in  Figure  7.  Because  of  this,  the  right  or  red  bar  will  have  the 
higher  amplitude  and  phase. 

The  screen  shot,  Figure  14,  was  taken  as  close  in  time  as  possible  to  the 
two  bars  at  the  horizontal  position.  Screen  shots  are  taken  while  the  program  is 
running  and  the  bars  are  in  motion. 
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Figure  14.  30  degree  incident  sound  with  bars  horizontal. 


Figure  15.  30  degree  incident  sound  with  bars  up  -  red  bar  higher  amplitude. 

It  is  difficult  to  see  in  these  transition  slides,  but  the  red  bar  also  has  phase 
difference  which  translates  to  a  slight  speed  difference  in  the  oscillations. 
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Figure  16.  30  degree  incident  sound  with  bars  down-  red  bar  higher  amplitude. 


3.  Sound  from  60  Degrees 

This  amplitude  and  phase  difference  is  even  more  noticeable  at  60 
degrees  incident  angle. 


Figure  17.  60  degree  incident  sound  with  bars  down-red  bar  highest  amplitude. 
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II.  BIOMIMETIC  MEMS  DESIGN 


Biomimetics  is  the  study  of  nature  and  the  attempt  to  mimic  that  nature  in 
man  made  devices.  After  all,  no  matter  what  religion  or  belief  system  one  might 
hold,  nature  has  at  least  a  few  thousand  year  head  start  on  the  production  of 
complicated  organisms  than  man  does.  This  time  advantage  has  most  often  lead 
to  incredibly  efficient  devices  worth  mimicry.  In  the  previous  chapter  it  was  shown 
that  biology  has  developed  a  sound  sensing  organ  capable  of  determining  the 
bearing  that  the  sound  came  from. 

MEMS  stands  for  Micro-Electro-Mechanical  System  and  “...  is  the 
integration  of  mechanical  elements,  sensors,  actuators,  and  electronics  on  a 
common  silicon  substrate  ...”  [MEMSNET,  2007],  Not  only  are  there  moving 
parts  at  the  micro  scale,  but  there  exists  the  ability  to  electrically  connect  to  those 
parts  so  that  a  sensor  can  be  built  and  readout  achieved.  Given  the  micrometer 
size  of  the  fly’s  ear  structure,  MEMS  is  the  perfect  medium  with  which  to  create 
the  biomimetic  sound  sensor. 

A.  MERITS  OF  THE  POLYMUMPS  PROCESS 

Previous  investigators  have  attempted  to  mimic  this  fly’s  hearing  organ  in 
silicon  devices  before;  for  examples  see  the  Yoo  et  al.  (2002),  and  Cui  (2004), 
papers.  However,  many  of  the  silicon-device  building  processes  were  custom 
processes  and  therefore  costs  reflected  those  specialized  processes. 

MUMPs®  stands  for  Multi-User  MEMS  Processes.  The  key  to  lowering 
the  cost  of  devices  is  to  try  to  use  one  of  the  existing  standard  processes. 
PolyMUMPs  is  one  of  those  existing  processes.  It  stands  for  Polysilicon  Multi- 
User  MEMS  Processes.  According  to  MEMSCAP,  the  owners  of  the  process, 
PolyMUMPs  is  a  “...  three-layer  polysilicon  surface  and  bulk  micromachining 
process,  with  two  sacrificial  layers  and  one  metal  layer...”  [MEMSCAP,  2003], 
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The  PolyMUMPs  process  helps  to  lower  the  cost  because  hundreds  of 
different  chips  can  be  integrated  in  the  same  batch  process.  In  addition,  the 
layers  of  metals  and  silicon  within  the  existing  process  definition  allow  for  easy 
electrical  access  to  the  moving  mechanical  parts  to  achieve  the  readout 
necessary  for  the  sensor  to  work. 

1.  PolyMUMPs  Process  Definition 

In  order  to  fully  understand  the  MEMS  device  design,  first  a  discussion 
about  each  of  the  different  layers  of  the  PolyMUMPs  process  is  required. 

Each  chip  starts  with  a  silicon  substrate.  The  substrate  acts  as  the 
foundation  on  which  the  structure  is  built.  It  is  approximately  400  micrometers  in 
depth. 

On  this  foundation  a  0.6  micrometer  of  nitride  is  grown.  The  purpose  of  the 
nitride  is  to  electrically  isolate  the  substrate  from  the  structure  above. 

Next  a  layer  of  0.5  micrometer  polysilicon  is  grown,  called  Poly  Zero. 
Polysilicon  is  silicon  in  which  there  is  no  preferred  orientation  in  which  the 
crystals  are  grown.  This  means  that  it  has  slightly  different  chemical 
characteristics  from  the  silicon  used  in  the  substrate  which  is  single  orientation 
crystalline  silicon.  Poly  Zero  layer  is  typically  used  for  electrical  contacts,  but  can 
be  designed  to  be  used  for  any  other  purpose.  In  the  case  of  the  sound  sensor 
built,  this  layer  was  used  for  electrical  contacts 

After  the  Poly  Zero,  a  layer  of  photo  resist  is  laid  down.  Photo  resist  is  a 
chemical  that  when  exposed  to  ultraviolet  light  will  change  its  chemical 
characteristics.  This  change  in  chemical  characteristics  allows  the  use  of  a 
pattern  mask  so  that  designs  can  be  placed  onto  the  photo  resist.  Then  the 
system  is  dipped  in  a  chemical  bath  or  reactive  ion  field  that  removes  the  photo 
resist  and  Poly  Zero  in  the  areas  not  exposed  to  light.  All  that  is  left  then  is  the 
patterned  Poly  Zero  and  the  exposed  photo  resist.  A  bath  in  a  different  solution 
removes  the  photo  resist. 
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Then  a  layer  of  Oxide  (Si02)  is  laid  down.  This  oxide  layer  is  two 
micrometers  thick.  The  oxide  is  patterned  using  photo  resist  and  a  mask 
following  the  same  basic  process  as  before.  The  purpose  of  the  oxide  layer  is 
usually  to  provide  a  temporary  space  and  later  will  be  chemically  removed  to 
allow  for  structure  to  be  free  released  and  able  to  move. 

The  next  step  is  another  layer  of  poly  silicon.  This  layer  is  called  Poly  One 
and  is  also  two  micrometers  thick.  Again  a  process  of  photo  resist,  masking,  and 
chemical  etching  a  design  takes  place.  In  the  built  sound  sensor  the  Poly  One 
layer  was  the  bottom  layer  of  the  bar-hinge-bar  system  mimicking  the  fly. 

On  top  of  the  Poly  One  layer  a  second  oxide  layer  is  placed.  Only  this 
oxide  layer  is  only  0.75  micrometers  thick.  It  follows  the  same  processes  as  the 
first  oxide  layer  and  is  typically  used  for  the  same  purposes. 

Following  patterning  of  oxide  two,  a  third  layer  of  polysilicon,  called  Poly 
Two  is  laid.  At  1.5  micrometers  thick  it  too  is  patterned  and  etched  following  the 
same  processes  described  before.  In  the  sound  sensor  this  layer  served  as  the 
top  layer  of  the  bar-hinge-bar  structure. 

The  final  layer  is  0.5  micrometers  of  gold  placed  on  top  of  the  Poly  Two 
layer  and  then  patterned  and  etched  in  similar  style  as  before.  The  gold  layer  is 
most  often  used  for  electrical  contacts. 

Figure  18  shows  the  layers  in  their  order  and  shows  thicknesses  as  well. 
The  substrate  and  nitride  are  not  shown,  as  they  are  only  the  foundation. 


Figure  18.  PolyMUMPs  layering  and  thicknesses  (After:  Koester  et  al.,  2003). 
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A  more  detailed  description  of  how  the  layers  in  the  PolyMUMPs  process 
are  created  can  be  found  in  the  PolyMUMPs  Design  Handbook  [Koester  et  al., 
2003],  However,  enough  of  the  basics  were  described  here  so  that  the  MEMS 
biomimetic  bar-hinge-bar  sound  sensing  device  can  be  discussed  in  detail  next. 

B.  BIOMIMETIC  MEMS  DESIGN 

In  the  fly,  the  tympanal  membrane  or  ear  drum,  translated  the  sound 
pressures  into  displacement  of  the  bars.  So  the  silicon  device  too  needed  to 
translate  from  sound  pressure  to  a  displacement.  The  schematic  of  the  MEMS 
structure  mimicking  the  fly’s  ears  is  shown  in  Figure  19. 


Figure  19.  Schematic  diagram  of  fly  ear  using  MEMS  (After:  Karunasiri  et  al., 
2005). 

The  basic  design  of  the  structure  using  PolyMUMPs  is  depicted  in  Figure 
20.  There  are  membranes  that  act  like  the  bars  and  the  tympanal  membranes  in 
the  fly’s  ears.  The  membranes,  which  are  often  referred  to  as  “wings”  since  they 
resemble  bird  wings  flapping  during  the  bending  mode,  are  attached  to  the 
substrate  via  a  support  block.  The  membranes  have  holes  in  them  attempting  to 
control  the  damping  of  the  system  as  the  wings  move  (Cs).  The  cantilever  beams 
connecting  the  support  block  and  the  wings  act  as  the  springs  Ks.  Again  Kt  and  Ct 
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are  from  the  material  properties  of  the  polysilicon  and  are  adjusted  by  the  amount 
of  silicon  placed  in  the  center  of  the  structure.  Initial  research,  by  Dr.  Kim,  proved 
that  a  design  of  1000  pm  by  2000  pm  size  plate  was  the  best  for  proper 
operation. 
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Figure  20.  PolyMUMPs  design  of  directional  sensor. 

What  is  difficult  to  see  in  Figure  20  is  the  small  air  gap  underneath  the 
device.  In  addition,  it  is  hard  to  see  what  material  each  part  is  made  from.  Figure 
21  shows  a  side  view  of  a  typical  device  showing  relative  depth. 
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Air  Gap  from  removing  Oxide  1  (2  micrometers) 


Nitride 

Substrate 

Poly  Zero 

(blue) 

(orange) 

Figure  21 .  Depth  view  of  a  typical  PolyMUMPs  design. 

When  the  wing  structure  moves,  the  small  air  gap  causes  a  phenomenon, 
called  squeeze  film  damping.  Squeeze  film  damping  is  basically  a  damping  force 
that  is  felt  on  the  wing  structure  due  to  the  air  trying  to  escape  from  underneath 
the  wing  as  it  moves  down.  If  no  holes  are  placed  in  the  plate,  then  this  damping 
force  can  be  very  strong  and  rendering  the  plate  immovable,  in  other  words,  the 
plate  does  not  react  to  sound.  If  there  are  too  large  of  holes,  or  too  many  holes, 
then  this  lack  of  damping  force  can  cause  the  structure  to  oscillate  too  much. 
Early  on  into  the  research,  it  became  apparent  that  controlling  the  damping  would 
be  important  to  causing  the  two  modes  of  operation,  rocking  and  bending,  to  be 
as  far  apart  and  distinct  as  possible.  Thus  multiple  device  designs,  each  with 
different  hole  sizes  and  number  of  holes  were  created.  Each  of  these  different 
designs  was  based  on  a  primary  design  with  slight  modifications  in  an  attempt  to 
find  the  optimum  configuration. 
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1. 


The  Primary  Design 


Figure  22  shows  the  primary  design.  The  device  is  smaller  than  the  device 
in  Figure  20,  because  additional  simulation  proved  that  for  a  3.5  micrometer  thick 
wing  structure,  a  device  approximately  1080  by  1080  micrometers  would  work 
better;  achieving  a  rocking  and  bending  frequencies  that  were  farther  apart,  with 
higher  amplitudes.  The  3.5  micrometers  comes  from  two  micrometers  of  Poly 
One  plus  1.5  micrometers  of  Poly  Two.  The  cantilever  beams  were  simulated  for 
a  length  of  75  micrometers  and  a  width  of  14  micrometers.  This  proved  to  be  a 
good  size  to  get  a  rocking  mode  of  approximately  3  khlz  and  a  bending  mode  of 
approximately  6.8  kHz.  The  holes  are  86  micrometers  square  and  there  are  four 
rows  of  nine  per  wing  structure.  The  small  holes  in  the  center  are  to  aid  in  Oxide 
layer  release. 


Figure  22.  Primary  PolyMUMPs  design. 

If  a  comparison  of  Figure  20  and  Figure  22  is  made,  it  is  quite  obvious  that 
a  few  additional  components  were  added  to  the  initial  design. 
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a.  Reasons  for  Additional  Components 


Figure  23  shows  the  electrodes  that  were  added  to  the  structure. 


They  are  made  out  of  the  Poly  Zero  layer.  The  purpose  of  the  electrodes  was  two 
fold:  First,  they  were  added  to  be  able  to  electrostatically  actuate  the  device  for 
proof  of  the  rocking  mode  and  bending  mode  operation  of  the  device.  This  was 
achieved  by  applying  an  electrical  signal  to  either  the  left  or  right  electrode,  and 
ground  to  the  center  electrode  which  is  connected  to  the  wings  via  Poly  One  and 
Poly  Two.  Results  of  this  testing  will  be  discussed  in  more  detail  in  Chapter  IV. 
Second,  although  not  originally  planned  for,  they  were  utilized  in  order  to  retrieve 
the  displacement  of  the  device  by  measuring  the  capacitance  of  the  wing-air- 
electrode  (plate-dielectric-plate)  of  the  right  or  left  side.  The  displacement  is 
proportional  to  the  capacitance.  A  more  detailed  discussion  of  this  capacitance 
measurement  system  will  be  discussed  in  Chapter  IV. 


Left  Electrical 
Connection  Pad 


Center 

Electrode 


Figure  23.  View  of  primary  device  electrodes. 
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Additionally,  there  were  diffraction  grating  “fingers”  added  to  the 
device;  a  set  of  moveable  grating  fingers  and  a  set  of  stationary  fingers.  Figure 
24  shows  the  gratings  up  close.  The  moveable  fingers  are  attached  to  the  left 
and  right  wings  and  will  move  with  the  wings  as  they  rock  and  bend.  The 
stationary  fingers  are  attached  to  their  own  support  structure  and  do  not  move. 


Figure  24.  View  of  grating  structure  attached  to  a  wing. 

The  purpose  of  the  diffraction  grating  is  to  measure  accurately  the 
displacement  of  the  wings.  A  near  infrared  laser  is  shined  on  the  fingers  and 
reflected  by  the  gold  on  top  of  the  grating  fingers.  The  diffraction  pattern  is 
measured  by  a  detector  as  shown  in  Figure  25.  The  displacement  can  be 
measured  down  to  less  than  a  nanometer  using  this  method. 


Laser 


Figure  25.  Workings  of  diffraction  grating  (From:  Kim  et  a!.,  2005). 
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b.  3D  View  of  Primary  Device  (without  Electrodes) 

Figure  26  shows  the  primary  device  in  three  dimensions.  This  view 
makes  it  easier  to  see  the  device.  The  electrodes,  which  are  normally  underneath 
the  wings  and  underneath  the  center,  were  removed  to  make  it  easier  to  see  the 
damping  holes.  The  red  is  the  polysilicon  layers  and  the  yellow  is  the  gold  layer. 


Figure  26.  3D  View  of  primary  device  without  electrodes. 

C.  CHIP  LAYOUT 

Proper  damping  in  this  device  is  very  crucial  to  the  success  of  accurately 
determining  the  bearing  of  the  sound  source.  Because  of  this,  it  was  necessary 
to  try  multiple  designs  to  physically  see  which  combination  of  hole  size,  number 
of  holes,  length  of  cantilever  beam,  etc.  provided  the  best  performance. 

Therefore,  21  different  designs  were  placed  on  a  single  chip.  Figure  27 
shows  the  layout  of  all  of  the  devices. 

Device  number  one  is  the  primary  design  discussed  earlier.  Deviating 
from  this  design,  device  two  has  larger  holes.  Device  number  three  has  smaller 
holes.  Testing  the  effect  of  beam  length  on  the  device,  device  number  four  has  a 
smaller  cantilever  beam.  Device  number  five  has  a  longer  cantilever  beam. 
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Combinations  of  changes  occurred  too,  for  example,  device  number  six  has 
larger  holes  and  a  shorter  beam.  These  changes  are  summed  up  in  Table  1 . 


Figure  27.  Chip  layout. 


DEVICE 

CANTILEVER  BEAM 

HOLES 

N  u  m  b  e  r 

Length 

(pm) 

Width 
(M  m  ) 

Length 

(Mm) 

W  id  th 
(Mm) 

Num  ber  per 
wing 

(across  x 
down) 

Square 

hole 

Length 

(Mm) 

i 

10  8  0 

10  8  0 

7  5 

1  4 

4x9 

8  6 

2 

10  8  0 

10  8  0 

7  5 

1  4 

3x7 

12  0  | 

3 

10  8  0 

10  8  0 

7  5 

1  4 

6x13 

5  0 

4 

10  8  0 

10  8  0 

3  7.5 

1  4 

4x9 

8  6 

5 

10  8  0 

10  8  0 

1  5  0 

1  4 

4x9 

8  6 

6 

10  8  0 

10  8  0 

3  7.5 

1  4 

3x7 

12  0  | 

i  7 

10  8  0 

10  8  0 

1  5  0 

1  4 

3x7 

1  2  0 

8 

10  8  0 

10  8  0 

3  7.5 

1  4 

6x13 

5  0 

9 

10  8  0 

10  8  0 

1  5  0 

1  4 

6x13 

5  0 

1  0 

8  5  0 

10  8  0 

7  5 

1  4 

3x9 

8  6 

1  1 

8  5  0 

10  8  0 

3  7.5 

1  4 

3x9 

8  6 

1  2 

8  5  0 

10  8  0 

1  5  0 

1  4 

3x9 

8  6 

1  3 

13  8  0 

10  8  0 

7  5 

1  4 

4x7 

1  2  0 

1  4 

7  8  0 

10  8  0 

7  5 

1  4 

2x7 

1  2  0 

1  5 

15  4  0 

10  8  0 

7  5 

1  4 

6x9 

8  6  | 

1  6 

15  5  0 

10  7  0 

7  5 

1  4 

9x13 

5  0 

1  7 

7  5  0 

10  7  0 

7  5 

1  4 

4x13 

5  0 

1  8 

10  8  0 

10  8  0 

7  5 

1  4 

4.5x9 
(holes  in 
center 
column  too) 

8  6 

1  9 

6  1  0 

6  1  0 

7  5 

1  4 

2x5 

8  6 

2  0 

10  8  0 

6  1  0 

7  5 

1  4 

4x5 

8  6 

2  1 

15  4  0 

6  1  0 

7  5 

1  4 

6x5 

8  6  | 

Table  1 .  Design  parameters  used  for  the  chip  in  Figure  27. 
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In  the  next  chapter,  the  simulation  of  the  devices  will  be  discussed  using 
the  previous  table  of  design  parameters. 
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III.  COMSOL  MODE  SIMULATION 


Simulations  of  the  sensors  designed  in  Chapter  II  were  conducted  both  in 
ANSYS  and  COMSOL  Multiphysics  finite  element  modeling  (FEM)  programs. 
ANSYS  simulations  were  conducted  prior  to  the  manufacturing  of  the  devices 
mostly  because  it  was  the  FEM  program  used  by  the  MEMS  department  at  Naval 
Postgraduate  School  at  that  time.  However,  during  the  project,  the  FEM  program 
COMSOL  Multiphysics  was  introduced  to  the  project  team.  COMSOL  is  much 
easier  to  use  than  ANSYS  and  provides  similar  results.  Therefore  all  simulations 
after  that  point  were  done  in  the  COMSOL  Multiphysics  program. 

The  ANSYS  simulations  were  completed  using  a  code  written  by  Dr. 
Byungki  Kim.  A  sample  of  this  code  is  contained  in  the  Appendix  B  for  reference. 

Given  the  comparable  results  between  the  two  FEM  programs,  and 
COMSOL’ s  ability  to  allow  the  user  to  easily  add  squeeze  film  damping  and 
sound  actuation  to  this  truly  multi-physics  problem,  made  COMSOL  the  clear 
choice  for  all  simulations. 

A.  COMSOL  MODE  SIMULATION  PROCESS 

All  the  COMSOL  Multiphysics  simulations  follow  the  same  basic 
processes.  The  process  is  as  follows:  draw  an  object,  assign  some  relevant 
physics  equations,  designate  the  materials  it  is  made  of,  establish  boundary 
conditions,  mesh  the  object,  set  the  solver  parameters,  and  then  finally  solve  the 
simulation.  In  many  cases  the  order  in  which  these  steps  are  conducted  is  not 
important,  but  it  makes  it  easier  not  to  skip  a  step  if  the  same  routine  is  followed 
every  time.  The  only  step  for  which  order  matters  is  the  solve  step,  it  must  be 
last. 

First,  either  a  2D  image  was  drawn  manually  or  imported  from  a  Computer 
Aided  Drafting  (CAD)  program.  Then  the  2D  image  was  extruded  to  a  3D  image 
with  the  proper  depth,  3.5  micrometers  in  this  case.  3.5  micrometers  was  the 
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sum  thickness  of  the  Poly  One  and  Poly  Two  layers,  which  comprised  the  wing 
structure.  Sometimes  drawing  in  2D  and  then  extruding  to  3D  is  much  easier 
than  attempting  to  draw  the  device  in  3D  from  the  start. 

After  extrusion,  the  relevant  physics  were  chosen.  In  this  case  since  the 
project  team  is  attempting  to  determine  the  Eigen  frequencies  of  the  structure, 
structural  mechanics  was  chosen;  specifically  the  MEMS  solid,  stress-strain 
option. 


After  selecting  the  physics,  the  construction  material  types  were  chosen. 
For  this  simulation  the  entire  structure  was  Polysilicon. 


Boundary  conditions  for  these  simulations  were  straightforward.  All 
portions  of  the  device  except  the  ends  of  the  cantilever  beams  were  free  to 
move.  The  ends  of  the  cantilever  beams  were  fixed,  illustrated  in  Figure  28. 


3  COMSOL  Multiphysics  -  Geom2/MEM5  Module  -  Solid,  Stress-Strain  (smsld) :  Casel  new.mph 


I  File  Edit  Options  Draw  Physics  Mesh  Solve  Postprocessing  Multiphysics  Help 


I  H  (?  a  A 1  Gpnmi  I  Geom2 1 

loundary  Settings  -  Solid,  Stress-Strain  (smsld) 

Fixed  ends  of 
Cantilever  Beams 


Figure  28.  Primary  device  with  fixed  boundary  conditions. 


After  setting  the  materials  and  boundary  conditions,  meshing  the  device 
was  important.  This  at  many  times  was  the  most  difficult  part.  If  it  is  meshed  too 
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finely,  then  the  solution  consumes  hundreds  of  megabytes  of  memory  and  the 
solution  time  increases  linearly  as  the  number  of  meshing  is  increased.  However, 
if  the  computer  system  does  not  have  enough  RAM,  the  COMSOL  program  will 
resort  to  using  virtual  memory  (hard  drive  storage)  and  in  this  case,  the  solution 
time  can  increase  by  one  to  three  orders  of  magnitude.  If  it  is  meshed  too 
coarsely,  then  the  physical  results  are  not  as  accurate.  Therefore  a  large  amount 
of  time  was  spent  trying  to  achieve  a  reasonable  meshing,  and  in  fact  a  solution 
that  was  mesh  independent.  In  the  case  of  simulating  just  the  Eigen  frequency 
modes,  this  was  not  as  important  an  issue  as  it  was  in  simulating  the  device 
actuation  under  a  sound  wave  which  will  be  discussed  in  Chapter  VI.  Figure  29 
shows  the  primary  device  mesh  with  the  normal  mesh  size  chosen  for  this 
simulation.  It  should  be  apparent  that  the  meshing  is  more  detailed  in  the  areas 
of  small  features  or  places  where  more  detailed  results  were  required. 


Figure  29.  Normally  meshed  primary  device. 
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Once  the  device  is  meshed  then  the  solver  parameters  can  be  set.  Since 
the  structure  in  Figure  29  has  two  Eigen  modes,  the  solver  parameters  were  set 
to  analysis  ->  Eigenfrequency,  solver  ->  Eigenfrequency,  and  the  desired 
number  of  Eigen  frequencies  were  set  to  three  (one  extra  for  exploring  higher 
order  modes),  as  illustrated  in  Figure  30. 


Figure  30.  Solver  parameters  settings. 

Finally  the  simulation  is  carried  out  to  determine  the  Eigen  frequencies. 
The  amount  of  time  that  the  simulation  takes  varied  greatly  depending  on  the 
device  complexity  that  was  simulated  and  the  tightness  of  the  meshing  chosen. 

1.  Rocking  and  Bending  Modes 

As  discussed  in  Chapter  I,  the  rocking  mode  looks  like  a  seesaw  as  it 
rocks  back  and  forth,  while  the  bending  mode  looks  like  a  bird  bending  its  wings 
or  flapping  its  wings  in  flight.  Figure  31  (a)  and  (b)  show  cartoon  time  sequenced 
representation  of  these  two  motions  using  the  simulated  data.  Figure  32  shows 
the  3D  view  of  the  primary  device  in  the  rocking  mode  and  Figure  33  shows  the 
primary  device  in  the  bending  mode. 
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Figure  31 .  Time  sequence  of  (a)  rocking  and  (b)  bending  modes. 
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Figure  32.  Primary  device  in  rocking  mode  representing  time  1  in  Figure  30  (a). 
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Figure  33.  Primary  device  in  bending  mode  representing  time  1  in  Figure  30 
(b). 

B.  COMPARISON  OF  ANSYS  AND  COMSOL  SIMULATION  RESULTS 

Table  2  shows  the  results  for  some  of  the  21  devices  simulated  in  ANSYS 
and  COMSOL  side  by  side. 
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DEVICE 

CANTILEVER 

BEAM 

HOLES 

ANSYS 

SIMULATED 

FREQUENCY* 

COMSOL 

SIMULATED 

FREQUENCY 

Number 

Length 

(Mm) 

Width 

(Mm) 

Length 

(Mm) 

Wdth 

(Mm) 

Square 

hole 

Length 

(Mm) 

Rocking 

Mode 

(kHz) 

Bending 

Mode 

(kHz) 

Rocking 

Mode 

(kHz) 

Bending 

Mode 

(kHz) 

1 

1080 

1080 

75 

14 

4x9 

86 

2.7 

6.8 

3.0 

6.8 

2 

1080 

1080 

75 

14 

3x7 

120 

2.9 

7.0 

3.1 

6.7 

3 

1080 

1080 

75 

14 

6x13 

50 

2.6 

6.8 

2.8 

7.0 

4 

1080 

1080 

37.5 

14 

4x9 

86 

3.7 

7.3 

3.8 

7.4 

5 

1080 

1080 

150 

14 

4x9 

86 

2.0 

5.8 

2.3 

5.7 

6 

1080 

1080 

37.5 

14 

3x7 

120 

3.8 

7.4 

3.8 

7.3 

7 

1080 

1080 

150 

14 

3x7 

120 

2.1 

6.0 

2.3 

5.7 

8 

1080 

1080 

37.5 

14 

6x13 

50 

3.8 

7.4 

3.5 

7.2 

9 

1080 

1080 

150 

14 

6x13 

50 

1.9 

5.6 

2.1 

5.8 

10 

850 

1080 

75 

14 

3x9 

86 

3.8 

7.9 

4.0 

7.7 

11 

850 

1080 

37.5 

14 

3x9 

86 

5.4 

8.8 

5.1 

8.6 

12 

850 

1080 

150 

14 

3x9 

86 

2.7 

6.5 

3.0 

6.3 

13 

1380 

1080 

75 

14 

4x7 

120 

2.1 

5.2 

2.3 

5.3 

14 

780 

1080 

75 

14 

2x7 

120 

4.2 

8.6 

4.5 

7.9 

15 

1540 

1080 

75 

14 

6x9 

86 

1.8 

4.5 

1.8 

4.7 

16 

1550 

1070 

75 

14 

9x13 

50 

1.7 

4.6 

1.7 

4.9 

17 

750 

1070 

75 

14 

4x13 

50 

4.2 

8.3 

4.3 

8.2 

18 

1080 

1080 

75 

14 

4.5x9 

86 

2.7 

6.2 

2.9 

6.4 

19 

610 

610 

75 

14 

2x5 

86 

6.1 

15.1 

6.3 

14.8 

20 

1080 

610 

75 

14 

4x5 

86 

3.3 

7.8 

3.5 

8.2 

21 

1540 

610 

75 

14 

6x5 

86 

2.1 

4.4 

2.2 

4.8 

*Simulated  with  Dr.  Kirris  ANSYS  code 


Table  2.  ANSYS  vs  COMSOL  simulations  for  the  chip  in  Figure  27. 


After  simulation,  testing  became  the  next  priority.  The  following  chapter 
discusses  the  details  of  the  lab  setup  and  the  testing  conducted. 
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IV.  DEVICE  TESTING 


The  directional  sensors  layout  using  MEMSCAP  software  (see  Figure  27) 
was  fabricated  by  the  MEMSCAP  foundry.  Figure  34  shows  an  actual  view  of 
one  of  the  sensors  from  an  optical  microscope.  The  grating  structures  at  the 
edges  of  the  wings  are  not  clearly  seen  due  to  the  low  magnification  of  the 
microscope.  However,  the  major  features  are  clearly  visible. 


Gratings 


Supporting  Block 


Right  Electrode 
Connection  Pad 


Left  Electrode 
Connection  Pad 


Supporting  Block 


Figure  34.  Optical  microscope  view  of  primary  device.  The  wing  structure  is 
1080  micrometers  by  1080  micrometers. 


A.  WIRING  AND  PACKAGING 

Using  a  Kyocera  68  Pin  Grid  Array  Package,  purchased  from  Addison 
Engineering  Inc.,  multiple  chips  were  packaged  and  wired  up  for  testing. 
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Wiring  was  conducted  using  a  4500  Digital  Series,  Kulicke  &  Soffa, 
Manual  Wire  Bonder  machine.  This  machine  solders  and  spools  wire  on  the 
micrometer  scale,  while  the  operator  controls  placement  via  a  microscope  and 
control  panel. 

Figure  35  shows  the  device  in  its  package  and  fully  wired.  In  Figure  35, 
the  chip  package  is  placed  in  a  retaining  hoop  to  mitigate  movement  of  the 
package;  this  was  part  of  the  testing  rig  that  will  be  explained  in  more  detail  later 
in  the  chapter. 


Figure  35.  Chip  in  its  package  and  wired. 

Figures  36  and  37  show  the  schematics  of  the  wiring  used  for  connecting 
a  set  of  sensors  located  in  the  outer  and  inner  parts  of  the  chip.  Figure  38  shows 
the  pin  numbering.  Pins  3,  8,  11,  and  15  were  used  to  supply  the  ground  signal 
for  sensors  in  Figure  36;  where  as,  pins  7  and  11  were  used  as  ground  for  the 
sensors  in  Figure  37.  Separating  the  chip  wiring  into  this  configuration  made  it 
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Outside  Chips 


easier  to  run  wires  as  each  wire  is  placed  manually  via  the  assistance  of  the 
machine.  If  all  21  designs  were  wired  on  one  chip,  then  it  would  have  been 
clustered  with  possibility  of  shorting  some  of  the  wires. 


Figure  36.  Outer  devices  wiring  diagram. 
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Inside  Chips 


Top 


Figure  37.  Inner  devices  wiring  diagram. 
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Figure  38.  Pin  numbering  diagram.  Pins  are  on  the  back.  Numbered  as  if  seen 
from  the  front.  Numbers  correspond  to  its  same  number  pad. 
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B.  EXPERIMENTAL  SETUP 


Figure  39  shows  a  picture  of  the  laboratory  setup  used  to  measure  the 
vibrational  amplitudes  of  the  sensors.  Each  of  the  21  devices  is  to  be  tested,  both 
electrostatically  and  via  sound  actuation,  in  order  to  test  for  their  performance. 


Figure  39.  Equipment  in  lab  setup  used  for  measuring  vibrational  amplitudes. 


1.  Functions  of  Equipment  Used 

As  Figure  39  shows,  there  were  some  key  components  integral  to  the 
testing  of  each  device.  They  were  the  reference  microphone,  the  sound  source, 
the  laser  vibrometer,  and  the  chip  being  tested.  Some  additional  equipment  was 
used  and  will  be  discussed  as  appropriate. 
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The  reference  microphone  was  placed  close  to  the  sensor  chip  to 
measure  the  sound  pressure.  It  is  a  Bruel  &  Kjaer  Pressure-field  1/8”  Microphone 
type  4138.  It  has  a  relatively  flat  response  curve  from  20  Hz  to  20  kHz  with  0.939 
mV/Pa  sensitivity. 

The  sound  source  was  a  Selenium  loudspeaker  type  DH200E  attached  to 
a  Hewlett  Packard  (HP)  331 4A  Function  Generator  and  HP  467A  Power 
Amplifier.  The  purpose  of  this  equipment  was  to  produce  a  set  amplitude  and 
frequency  of  the  incident  sound  wave.  Incident  sound  wave  angle  was  adjusted 
manually  by  adjusting  the  speaker  position. 

The  laser  vibrometer  was  a  Polytec  single  point  vibrometer  model  OFV 
302,  with  a  model  OFV  2600  controller.  Its  purpose  was  to  measure  the 
displacement  of  the  wings,  to  a  resolution  in  the  sub  nanometer  range. 

Additional  equipment  not  shown  in  Figure  39,  but  just  as  valuable  was  a 
lock-in  amplifier  and  an  oscilloscope.  The  lock-in  amplifier  was  and  EG&G 
Princeton  Applied  Research  model  5210.  Lock-in  amplifiers  are  used  to  measure 
the  amplitude  and  phase  of  repetitive  AC  signals  buried  in  noise.  The 
oscilloscope  was  a  four  channel  Agilent  Infiniium  DSO8064A  model.  It  was  used 
to  show  the  signal  from  the  reference  microphone  and  the  device  being  tested  on 
the  same  screen. 

C.  ELECTROSTATIC  TESTING 

Prior  to  Naval  Postgraduate  School  receiving  their  own  laser  vibrometer, 
Teledyne  Technologies,  Inc.  (formerly  Rockwell  Scientific)  in  Thousand  Oaks, 
CA,  allowed  the  project  team  to  use  their  laser  vibrometer. 

Testing  was  conducted  at  their  facility  in  Thousand  Oaks,  and  therefore 
needed  to  be  limited  in  the  number  of  devices  tested.  The  devices  chosen  to  test 
were  4,  5,  8,  14  and  19  as  indicated  in  Figure  37.  Each  device  was 
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electrostatically  actuated  at  a  voltage  of  one  volt  over  a  frequency  sweep  from 
zero  to  20  kHz.  Note  that  for  actuating  the  rocking  mode  the  bias  should  be 
applied  to  only  one  side  of  the  membrane. 

Testing  started  with  device  4  because  it  was  easy  to  access  and  easiest  to 
manually  wire  the  pins  to  the  electrostatic  source.  In  addition  device  4  was  similar 
to  the  primary  device  except  it  had  a  smaller  cantilever  beam  size.  Results  from 
device  4  did  not  show  two  resonant  frequency  amplitude  peaks  as  were 
expected.  Results  are  shown  in  Figure  40;  note  that  the  Y  axis  is  velocity 
amplitude.  Remember,  the  bending  mode  was  predicted  to  be  at  approximately 
7.5  kHz  for  this  device,  which  is  close  to  the  measured  value  in  Figure  40. 
However,  no  rocking  mode  was  seen  on  the  graph.  The  peak  near  20  kHz  is  a 
higher  order  mode. 


Figure  40.  Device  4,  magnitude  of  velocity  versus  frequency  curve. 

This  was  perplexing,  and  because  of  this  result,  device  8  was  selected  to  be 
tested  next.  It  was  assumed  that  perhaps  the  missing  rocking  mode  was  due  to 
the  damping  which  broadens  the  peak  reducing  the  amplitude. 
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Figure  41 .  Device  8,  magnitude  of  velocity  versus  frequency  curve. 

Device  8  was  the  same  dimensions  as  device  4,  in  length  and  width  of  the 
wings,  and  also  in  cantilever  beam  size;  the  differences  were  the  hole  size  and 
number  of  holes.  Device  8  did  not  show  the  desired  rocking  mode  either,  as 
shown  in  Figure  41 . 

Continuing  to  be  perplexed,  the  affects  of  the  cantilever  beam  length  were 
probed  next  by  selecting  device  number  5  as  the  next  device  to  be  tested.  Again, 
results  were  not  as  expected.  The  longer  cantilever  beam  did  not  help  the 
rocking  mode  to  appear. 

At  this  point,  it  was  decided  to  test  a  device  with  smaller  wing  dimensions. 
Remember  that  the  rocking  mode  frequency  was  inversely  proportional  to  mass 
(equation  1 .20).  The  prevailing  thought  was  that  a  smaller  device  would  have  a 
higher  rocking  frequency  and  therefore  might  allow  the  rocking  mode  to  be  seen 
more  clearly.  Following  this  train  of  thought  device  14  was  chosen  next.  Device 
14  did  not  show  the  desired  results  either,  but  device  14  was  only  smaller  in 
length. 

Finally,  device  19  which  had  the  smallest  wing  dimensions  was  tested. 
Device  19  provided  the  frequency  response  curve  that  was  expected.  Figure  42 
clearly  shows  the  rocking  and  bending  mode  amplitude  peaks  as  expected. 
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Figure  42.  Device  19,  magnitude  of  velocity  versus  frequency  curve. 


If  the  results  are  converted  to  displacement  amplitude  vice  velocity 
amplitude  (displacement  =  velocity/ co)  ,  the  results  are  4.75  kHz  and  15.1  kHz. 
Remember,  from  Table  2,  the  expected  rocking  mode  was  approximately  at  6.0 
kHz  and  the  bending  mode  was  expected  at  approximately  15  kHz.  However, 
these  results  of  approximately  4.75  kHz  and  15.1  kHz  were  within  reason  and 
promising.  Explanation  as  to  the  differences  between  simulation  and  testing  is 
discussed  in  Chapter  V. 


1.  Visual  Results 


Besides  the  frequency  response  curves  shown  earlier,  the  laser 
vibrometer  computer  program  that  accompanied  the  controller  allowed  the 
operator  to  superimpose  an  optical  image  over  an  animation  of  the  actual  device 
amplitude  data.  An  animation  can  easily  be  created  during  this  process.  Figure 
43  shows  the  rocking  mode  of  device  19  at  4.75  kHz  in  this  manner.  Figure  43  is 
shown  in  three  screen  shots  side  by  side,  as  the  device  rocks  from  the  left  side 
down  to  the  right  side  down.  Similarly,  Figure  44  shows  the  bending  mode  at 
15.1  kHz  as  well.  These  figures  are  a  graphical  representation  of  the  actual  data. 
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Figure  43.  Device  19  electrostatically  driven  rocking  mode. 


Figure  45  shows  a  detailed  image  of  device  19  at  one  instance  of  time 
while  it  is  undergoing  a  rocking  mode.  The  number  of  data  points  taken  by  the 
laser  on  the  gratings  was  minimal  to  save  in  computation  time,  thus  some  noise 
is  seen  in  the  figure;  however,  the  extreme  detail  was  useful  in  some  instances. 


52 


Figure  45.  Device  19  electrostatically  driven,  at  4.75  kHz,  rocking  mode 
detailed  view. 


2.  Computer  Controlled  Evaluation  Board  with  M3110 
Capacitance  Measurement  Chip 

Early  on  in  the  electrostatic  testing  phase,  it  became  apparent  the  Dr. 
Kim’s  optical  method  of  measuring  displacement  would  not  work  well  with  these 
PolyMUMPs  devices.  It  turns  out  that  due  to  the  length  of  the  gratings  and  the 
manufacturing  process  that  there  is  residual  stress  in  the  gratings  that  cause  the 
gratings  to  bend  slightly  upward.  The  net  result  of  this  bend  is  that  the  optical 
laser  system  discussed  in  Chapter  II  and  shown  in  Figure  27,  did  not  allow  the 
diffraction  pattern  to  accurately  determine  displacement.  A  different  system 
needed  to  be  employed. 

Many  MEMS  devices  in  industry  use  the  capacitive  approach  to 
measuring  displacement  and  a  commercial  board  was  readily  available  to 
achieve  the  accuracy  required  for  this  device.  The  stated  sensitivity  of  the  board 
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is  approximately  four  aF/(Hz)1/2,  remarkably  sensitive  [www.irvine-sensors.com]. 
In  addition  each  of  the  techniques  used  on  this  board  could  eventually  be 
miniaturized  and  easily  placed  on  the  same  chip  as  the  biomimetic  sound  device 
being  tested. 

Figure  46  shows  the  Irvine  Sensors,  MS3110BDPC  board  attached  to  a 
chip  about  to  be  tested  electrostatically;  the  electrostatic  connections  have  not 
yet  been  made  though. 


Figure  46.  Computer  controlled  evaluation  board  with  M3110  capacitance 
measurement  chip. 

D.  SOUND  TESTING 


Once  the  concept  of  rocking  modes  and  bending  modes  of  the  actual 
biomimetic  devices  were  proven  with  electrostatics,  it  was  imperative  to  test  the 
devices  in  a  sound  field. 
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1.  Polytec  Laser  Vibrometer 


The  Naval  Postgraduate  School  received  its  Polytec  single  point 
vibrometer  model  OFV  302,  with  a  model  OFV  2600  controller  and  the  lab  was 
setup  as  shown  previously  in  Figure  39.  The  focus  was  on  device  19  since  this 
device  showed  both  a  clear  rocking  and  bending  mode  during  electrostatic 
testing  and  the  most  data  was  available  for  validation  of  the  results.  The  sound 
pressure  at  the  device  was  about  1  Pascal  (approximately  95  dB). 

2.  Results 


Figure  47  shows  the  response  of  device  19  when  actuated  using  a  sound 
wave.  It  is  showing  the  bending  mode  at  15.1  kHz;  again  it  is  three  screenshots 
side  by  side  to  portray  a  time  sequence.  Displacement  was  found  to  be  less  than 
one  nanometer  primarily  due  to  leaking  of  the  sound  pressure  through  the  holes 
used  to  control  the  damping.  The  rocking  mode  amplitude  was  too  small  to  be 
measured  since  it  depends  on  the  pressure  difference  between  the  two  wings.  A 
detailed  explanation  is  given  in  Chapter  V  explaining  this  behavior. 


Figure  47.  Actuation  of  device  19  using  sound  incident  normal  to  the 
membrane  surface.  The  bending  mode  can  be  clearly  seen  from  the 
time  series. 
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E.  ADDITIONAL  INFORMATION  OBTAINED  FROM  FREQUENCY 

RESPONSE  CURVES 

Figure  48  shows  the  frequency  response  curve  from  one  of  the  testing 
runs  on  device  19.  It  clearly  shows  a  rocking  mode  peak  and  a  bending  mode 
peak.  Additional  information,  though,  can  be  gained  from  such  graphs.  The  full- 
width-at-half  maximum  of  the  resonance  peak  (normal  mode)  can  be  used  to 
determine  the  quality  factor  (Q).  Q  is  defined  as  the  ratio  of  resonance  curve 
width  at  half  maximum  to  the  peak  resonance  frequency  ( co0 ).  The  full-width-at- 

half  maximum  (FWHM)  is  approximately  equal  to  two  times  the  damping  constant 
(/?)  [Taylor,  2005],  That  means  the  Q  value  is  a  qualitative  measurement  and  is 
inversely  proportional  to  the  damping  constant,  as  seen  in  equation  1 .34  [Kinsler 
etal.,  2000], 

Q  =  coJ  2j8  (1.34) 

Therefore  by  looking  at  the  frequency  response  curves,  a  clear  indication 
of  the  damping  can  be  seen.  The  wider  the  resonance  curves  at  the  half 
maximum  points,  the  stronger  the  damping.  Additionally,  the  damping  coefficient 
can  be  calculated  within  some  minor  margin  of  error,  directly  off  the  response 
curves  and  thus  evaluated  device  changes  more  rapidly. 

Equation  1 .34  assumes  that  the  amplitude  of  the  resonance  peak  is  power 
amplitude;  however,  the  laser  vibrometer  output  is  velocity  amplitude  as  shown  in 
Figure  48.  Therefore,  a  formula  that  utilized  the  velocity  amplitude  was  derived. 
Since  full-width-at-half  maximum  (FWHM)  applies  to  power  which  is  proportional 
to  amplitude  squared,  full-width-square  root  of  2-maximum  [FWJlM]  applies 
here.  Derivation  is  shown  in  Appendix  C,  but  the  final  result  is  equation  1.35. 

Fw42M*2/3  (1.35) 

Sometimes  it  is  easier  to  measure  the  full-width  at  half-maximum,  and 
therefore  using  a  similar  derivation  it  can  be  shown  that  FWHM  =  V3  *2/? . 
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The  damping  parameter  can  be  estimated  using  either  approach.  Then  a 
comparison  of  the  damping  constant  to  that  of  the  modal  frequency  peak  is 
compared  using  the  relationships  listed  as  equation  1.36  [Taylor,  2005], 

(3<co0 =>  under  damped 

{3  =  co0^>  critically  damped  (1 .36) 

ft  >  coQ  -=>  over  damped 

A  comparison  of  the  damping  constant  to  that  of  the  rocking  mode 
frequency  for  device  19  using  the  data  shown  in  Figure  48  was  conducted  and 
the  result  is  that  device  19,  with  its  large  holes,  is  slightly  under  damped; 

179691^ <270 17  —  ^  fd<C00^  under  damped.  Keep  in  mind  that  there 

sec  sec 

is  only  a  two  micrometer  air  gap  under  the  wing  structure  and  most  of  that  air 
needs  to  move  through  the  relatively  large  holes  in  order  for  an  appreciable 
displacement  to  occur,  thus  slightly  under  damped  makes  sense. 


Figure  48.  Velocity  as  a  function  of  excitation  frequency  for  device  19. 
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V.  POSSIBLE  REASONS  FOR  DIFFERENCES  BETWEEN 
SIMULATIONS  AND  EXPERIMENTAL  RESULTS 


The  first  question  that  comes  to  mind  after  seeing  the  testing  results  is 
why  did  the  simulation  predict  a  rocking  mode  at  approximately  6  kHz  and  the 
testing  result  in  a  rocking  mode  of  4.75  kHz?  The  bending  mode  seemed  to  be 
almost  right  on,  as  the  simulation  predicted  approximately  15  kHz  and  the  testing 
showed  15.1  kHz.  What  could  cause  the  differences? 

A.  TOLERANCES  IN  STRUCTURE  SIZE 

The  first  insight  came  while  reviewing  the  optical  microscope  images  more 
closely.  It  seems  that  the  tolerances  in  the  fabrication  process  allowed  for  minute 
changes  to  the  dimensions  of  the  devices. 


Figure  49  shows  an  optical  microscope  view  of  device  19  with 
measurements  of  the  some  of  the  key  dimensions  that  control  the  rocking 
frequency. 


Figure  49.  Device  19  optical  microscope  view  with  measurements. 
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Remember,  from  Table  2,  the  torsion  beam  was  simulated  to  be  75 
micrometers  long  and  14  micrometers  wide.  Figure  49  clearly  shows  that  the 
cantilever  beam  on  the  actual  device  tested  was  0.4  micrometers  longer  and  0.1 
shorter  in  width. 


B.  SIMULATION  SIMPLIFICATIONS 

In  addition  to  variation  of  the  torsion  beam  dimensions,  the  gold  layer  on 
the  polysilicon  was  not  incorporated  in  the  earlier  simulation.  Recall,  equations 
1.20  and  1.21  showed  that  the  modal  frequencies  were  inversely  proportional  to 
the  mass  of  the  device.  Not  including  the  gold  might  significantly  affect  the  mass. 
So  in  order  to  incorporate  the  effects  of  the  gold  layer,  additional  simulations  in 
COMSOL  were  conducted.  The  result  is  the  rocking  mode  was  lowered  to  4.9 
kHz  which  is  much  closer  to  the  4.75  kHz  achieved  during  testing. 

1.  Revised  COMSOL  Simulation  Results 

Figure  50  shows  the  rocking  mode  modal  frequency  at  4.9  kHz.  If  Figure 
49  is  compared  to  Figure  50,  the  gold  strips  stand  out  more. 

eigfreq_smsld(  1)=4918.43 1044  Boundary:  Total  displacement  [nm]  Deformation:  Displacement  [m] 

< 


a- 


Figure  50.  Device  19  COMSOL  simulation  with  gold  layer. 

60 


C.  HOLE  SIZE  TOO  LARGE  CAUSING  SOUND  LEAKAGE 

As  was  mentioned  in  Chapter  IV,  displacement  of  device  19  when 
subjected  to  an  incident  sound  field  was  minimal  at  less  than  one  nanometer. 
This  was  perplexing.  Simulations  had  shown  that  displacement  at  the  rocking 
modal  frequency  should  have  been  in  the  tens  of  nanometers  at  the  smallest. 
True  these  simulations  were  not  simulations  of  the  device  subjected  to  a  sound 
field,  but  it  was  believed  that  if  the  modal  frequency  was  reached,  displacement 
amplitudes  should  have  been  comparable.  So  why  is  it  that  the  displacement 
was  so  low? 

During  the  design  phase,  it  was  believed  that  as  long  as  the  damping 
holes  were  smaller  than  the  wavelength  of  the  driving  sound  that  the  holes  would 
be  negligible  and  the  device  would  operate  similar  to  that  of  a  radar  dome.  In  the 
radar  dome,  the  holes,  do  not  affect  the  electromagnetic  properties  of  the  dome. 

The  key  physics  mistake  made  here  is  that  electromagnetic  waves  are 
transverse  which  actively  interacts  with  the  meshed  surface  of  the  radar  system; 
whereas  sound  waves  are  longitudinal  and  pass  through  the  holes  which  reduce 
the  pressure  difference  required  for  driving  the  wing  membrane.  Therefore  with 
too  large  of  holes  in  the  membrane,  it  is  very  easy  for  the  pressures  on  either 
side  of  the  membrane  to  nearly  equalize;  thus  minimizing  displacement. 

1.  Solid  Plate  Membrane 

Simulations  were  conducted  for  a  device  with  wing  sizes  comparable  to 
that  of  device  19,  but  without  the  holes.  Device  19  is  610  micrometers  square 
which  equates  to  372,100  square-micrometers  in  area. 

This  solid  plate  device  was  simulated  in  a  sound  field  of  amplitude  1  Pa,  at 
approximately  8.5  kHz.  Displacement  was  approximately  70  nanometers  for  the 
rocking  mode  as  shown  in  Figure  51 . 
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sound  actu-ated  45deg  left  to  right 
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Figure  51.  Solid  Plate,  sound  actuated  45  degrees  left  to  right,  COMSOL 
simulation  screen  shot. 

2.  Affect  of  Holes  in  the  Membrane 

Now  that  it  was  known  that  a  solid  plate  achieves  a  reasonable 
displacement,  it  was  decided  to  study  exactly  what  affect  holes  had  on  the  overall 
displacement  of  the  device  19.  Figure  52  shows  the  displacement  of  the  wings 
for  device  1 9  under  the  same  sounds  field  used  for  simulation  of  the  solid  plate  in 
Figure  51.  The  maximum  displacement  in  this  case  was  found  to  be  less  than  a 
nanometer  compared  to  the  almost  70  nanometers  for  the  solid  structure.  The 
measured  velocity  (v)  data  for  device  19  gave  about  8  micrometers/sec  at  the 
peak  of  the  resonance  for  rocking  mode.  Since  the  velocity  and  displacement  (x) 
amplitudes  are  related  by  v  =  cox,  the  estimated  displacement  amplitude  using 

the  measured  velocity  data  for  the  rocking  mode  is  about  0.3  nanometers.  This  is 
the  same  order  of  magnitude  as  the  simulated  value  considering  the  fact  that  the 
simulation  gives  an  over  estimation  due  to  the  exclusion  of  the  damping  effect. 
The  Table  3  summarizes  the  simulated  and  experimental  displacement 
amplitudes  for  device  19  along  with  that  of  the  solid  structure.  The  reason  for  the 
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difference  in  the  resonance  frequencies  is  that  the  simulated  device  does  not 
contain  the  grating  fingers,  as  the  do  not  affect  the  displacement,  and  only  added 
additional  calculation  time  and  memory  to  the  simulation. 


Bou  ndary:  T  otal  di  splacem  ent  [m  ]  D  ef  orm  ation :  D  ispl  acem  ent  [n 


Max:  9.608e-10 


Min:  0 


Figure  52.  COMSOL  simulation  screen  shot  of  the  device  19  under  the  same 
sound  field  as  the  solid  plate  in  Figure  51 . 


Resonance  Frequency 

Area  Removed 

Device  Displacement 

Solid  Plate 

372,100(pm2) 

610  pm  x  610  pm 

8,490  Hz 

none 

67  nm 

Simulated  device  19 

(without  the  grating 

fingers) 

11,167  Hz  (without  the 

grating,  lower  mass  =>  higher 

resonance  frequency) 

147,920  pm2  (20  holes 

each  86  pm  x  86  pm) 

0.9  nm 

Measured  device  19 

4.75  Hz 

Same  as  above 

0.3  nm 

Table  3.  Comparison  of  simulated  and  experimental  displacements. 
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D. 


INCORPORATION  OF  SQUEEZE  FILM  DAMPING  INTO  THE 
SIMULATIONS 


Squeeze  film  damping  is  a  type  of  damping  that  occurs  when  a  thin 
volume  of  a  gas  exists  under  a  moving  plate.  As  the  plate  moves  down  it  needs 
to  displace  the  air  underneath;  damping  comes  from  the  rate  at  which  the  plat 
can  displace  this  gas.  The  damping  reduces  the  vibration  amplitude  and 
broadens  the  resonance  peak.  This  is  the  primary  reason  the  simulation  did  not 
match  the  experimental  data  since  the  squeeze-film  damping  was  not  included  in 
early  simulations.  Around  the  time  that  the  discrepancy  between  the  simulation 
and  experimental  data  was  discovered,  COMSOL  released  its  version  3.3a 
version.  Built  into  this  version  was  the  ability  to  simulate  squeeze-film  damping 
in  a  way  that  did  not  appreciably  increase  the  memory  needs  or  the  simulation 
time. 


1 .  COMSOL  Process 

Section  7,  of  the  COMSOL  “MEMS  MODULE”  manual,  discusses  the 
details  of  the  theory,  but  it  is  based  on  a  modified  Reynolds  equation  [COMSOL, 
2006],  This  whole  process,  although  not  COMSOL  specific,  is  also  described  in 
the  Veijola  et  al.  2005  paper.  Figure  53  (a)  shows  visually  that  the  air  is  being 
forced  out  from  under  the  plate  and  escapes  through  the  perforated  holes  in  the 
plate.  Figure  54  shows  the  equations  that  COMSOL  is  using  in  its  film  damping 
physics  boundary  condition  section;  Figure  54  also  shows  many  of  the  settings 
used  for  simulations  for  the  regions  defined  as  holes.  The  regions  defined  as  not 
holes,  have  its  relative  diffusivity  (Dh)  set  to  1,  its  relative  compressibility  (Ch)  and 
its  Perforation  admittance  (Yh)  set  to  zero  as  shown  in  Figure  55.  Figure  53  (b) 
shows  graphically  what  the  settings  should  be  and  why. 
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Figure  53.  Film  damping  settings  based  on  Veiiola  et  at.  2005  paper  (From: 
Veijola  et  ai,  2005). 


Figure  54.  Film  damping  settings  for  the  regions  defined  as  holes. 
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Figure  55.  Film  damping  settings  for  the  regions  not  defined  as  holes. 

2.  Device  #19  Results 

Figure  56  shows  that  device  19  is  under  damped  as  was  calculated 
previously  in  Chapter  IV.  The  simulation  was  to  apply  a  force  of  one  Pascal  to 
the  left  wing  (which  is  the  upper  half  in  the  figure),  and  determine  the  damping 
pressure  magnitude.  Figure  56  shows  that  there  is  an  average  damping 
pressure  less  than  0.5  Pascals.  Even  though  there  are  a  few  places  with  spots  as 
high  as  1.5  Pascals  the  average  damping  is  less  than  0.5  Pascals.  Further 
research  is  needed  to  validate  these  findings  and  to  ensure  that  the  COMSOL 
simulation,  using  the  modified  Reynolds  equation,  is  accurate.  However,  based 
on  the  frequency  response  curve  damping  calculations  conducted  in  Chapter  IV, 
which  showed  an  under  damped  situation,  the  results  of  the  simulation  are  in 
agreement  that  the  device  is  under  damped. 
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Figure  56. 


Device  19  damping  pressure  results  for  actuation  at  6.1kHz  for  a 
one  Pascal  force. 
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VI.  CONCLUSION 


It  was  shown  that  it  is  possible  to  design  and  fabricate  a  biomimetic 
MEMS  based  directional  acoustic  sensor  that  displays  the  expected  rocking  and 
bending  modes  of  vibration,  similar  to  those  of  the  fly  Ormia  ochracea. 
Fabrication  of  the  sensor,  utilizing  the  PolyMUMPs  process  was  also  proven 
possible.  However,  the  current  designs  evaluated  produced  amplitudes  of 
vibrations  one  order  of  magnitude  smaller  than  expected  values.  In  which  case, 
simulations  agreed  well  with  obtained  results.  Thorough  simulation  studies 
indicated  that  due  to  the  longitudinal  character  of  sound  waves,  the  front  side 
holes  used  for  reducing  damping  greatly  diminished  the  pressure  difference 
across  the  plate  surfaces.  In  which  case,  simulations  agreed  well  with  obtained 
results.  Further  efforts  were  performed  to  design  a  working  model.  Simulations 
showed  that  a  solid  plate  would  produce  a  larger  displacement  in  the  sound  field. 
Simulations  showed  that  a  solid  plate  would  produce  a  larger  displacement  in  the 
sound  field.  However,  damping  becomes  a  large  issue  with  a  solid  plate.  Further 
modeling  efforts  showed  encouraging  results  where  by  locating  holes  in  the 
substrate  (rather  than  in  the  vibrating  plates)  the  damping  can  be  controlled  to 
desired  levels.  However,  this  direction  leads  to  specialized  processes  -  thus 
increasing  costs.  Nonetheless  recent  modeling  efforts  on  Silicon-On-Insulator 
(SOI)  MUMPs  process,  which  allows  for  etching  of  the  substrate,  indicate 
potential  for  designing  this  biomimetic  MEMS  device  on  a  standardized  process. 
Further  work  is  required  to  obtain  an  optimum  design  for  this  biomimetic  MEMS 
sensor  that  can  be  fabricated  on  a  standard  process. 
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VII.  RECOMMENDATIONS  FOR  FUTURE  WORK 


A.  INVESTIGATION  OF  SOIMUMPS  PROCESS 


In  order  to  achieve  the  optimized  performance  of  the  sensor  it  is 
recommended  that  future  research  be  focused  on  the  SOIMUMPs  process, 
focusing  on  solid  plates.  Damping  control  can  be  achieved  via  batch  etching  the 
damping  holes  in  the  substrate  for  minimal  cost  increase.  Figure  57  shows  a 
graphical  idea  of  what  the  sensor  should  look  like  in  3D;  the  plate  is 
approximately  10  micrometers  thick,  while  the  substrate  is  400-500  micrometers 
thick.  Figure  58  shows  the  same  SOI  device,  but  is  a  side  view  so  that  the 
etched  holes  are  visible. 


400-500  Micrometers 
Thick  Substrate 


Metal  for  Capacitive 

Displacement 

Measurements 


Solid  Plate 


Damping 

Holes 


Figure  57.  Recommend  Sensor  layout  with  perforated  holes  on  the  back  plate 
to  increase  vibration  amplitude  and  reduce  squeezed  film  damping. 
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Wings  with  electrodes 


Air  gap 


Perforated  holes 


Figure  58.  Recommend  Sensor  layout  side  view  showing  etched  holes  under 
wing  structure.  Again  substrate  is  400-500  micrometers  thick,  while 
the  wings  are  only  10  micrometers  thick,  so  this  is  not  drawn  to 
scale. 


B.  COMSOL  SIMULATIONS  OF  THE  SOIMUMPS  PROCESS 

Finally  simulations  should  be  done  for  how  the  air  will  escape  through  the 
etched  holes  in  the  substrate  and  the  damping  that  is  created  from  this  air  flow. 

C.  ADDITIONAL  WORK 

There  is  much  room  for  future  research  related  to  this  Thesis  topic.  Listed 
below  are  just  a  few  examples  of  future  research  needed  to  perfect  this 
biomimetic  sensor,  not  related  to  the  SOIMUMPs  process.  They  are  in  no 
particular  order. 

a)  Verify  the  COMSOL  simulation,  using  the  modified  Reynolds 
equation,  is  accurate. 

b)  Working  out  the  optimum  settings  for  the  M3110  Capacitance 
Measurement  Chip  and  the  placement  of  the  boards  components 
on  the  same  chip  as  the  sensor. 

c)  Determining  power  consumption  and  bearing  resolution. 

d)  Verify  operation  in  a  noisy  environment;  measure  signal  to  noise 
ratio.  Determine  the  degradation  in  the  bearing  resolution  in  this 
noisy  environment. 

e)  Working  out  placement  of  devices  to  measure  in  2D,  and  then 
eventually  in  3D. 
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f)  Packaging  issues  for  device  to  work  underwater;  keeping  in  mind 
the  impedance  differences  between  the  water  and  the  air  in  the 
package  might  affect  the  angle  information. 

g)  Creating  an  array  of  perfected  sensors. 
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APPENDIX  A  -  MATLAB  CODE  SAMPLE  FOR  FLY  ANIMATION 


o, 

o 

o, 

o 

o, 

o 

Q, 

O 

O, 

O 

g. 

o 

o, 

o 

o, 

o 

o, 

o 

o, 

o 

Q. 

o 

o, 

o 

g, 

o 

o, 

o 

o, 

o 

g. 

o 

o, 

o 

g„ 

o 

g. 

o 

o, 

o 

g. 

o 

o, 

o 

g„ 

o 

g. 

o 

o, 

o 

g. 

o 

o, 

o 

g„ 

o 

g. 

o 

o, 

o 

g. 

o 

o, 

o 


Script  file:Fly_ear  animatedrev2 .m 

Purpose:  Program  for  calculating  displacement  of  the  two  membranes 
Values  are  referenced  from  the  1995  Miles  et  al:  Mechanically  coupled 
ears ....  paper . 

The  program  does  the  following: 

1 .  Calculates  the  membrane  displacements  left  and  right 

2 .  Animates  the  results 

Personnel  who  have  worked  on  the  program: 

1)  LT  Tim  Shivok 

2)  Jeff  Catterlin  (high  school  intern) 

3)  Professor  Guest 


Record  of 

Date 

revisions : 

Programmer 

Class 

Description  of  change 

08/15/06 

Timothy  Shivok 

Thesis 

Original  code 

08/31/06 

Prof.  Guest 

revl-Modif ied  code  to  add 
animated  results 

09/05/06 

Timothy  Shivok 

Thesis 

rev2-Modif ied  code  to 
solve  for  wt,wr,et,  and  er 
vice  use  set  results.  Also 
added  the  calculations  for 
both  X  and  Y  values  vice 
just  one. 

09/07/06 

Prof.  Guest 

Added  code  to  create  avi 
movie 

clear  all; 
clc; 

o, 

o 

%  Set  the  constants 

m=2.88e-10;  %kg  -->  effective  mass  of  the  organ 

kl=.576;  %  N/m  -->  left  membrane  spring  constant 

k2=.576;  %  N/m  -->  right  membrane  spring  constant 

k3=5.18;  %  N/m  -->  coupling  of  membrane  spring  constant 

cl=1.15e-5;  %N  s/m  -->  left  membrane  damping  constant 

c2=1.15e-5;  %N  s/m  -->  right  membrane  damping  constant 

c3=2.88e-5;  %N  s/m  -->  coupling  of  membrane  damping  constant 

s=2.88e-7;  %m^2  -->  surface  area???????? 

d=1.2e-3;  %m  -->  distance  between  the  force  locations 

g. 

o 

%  Calculate  freq  and  damp  ratio 

g. 

o 


wr= (kl/m) A . 5;  %calculates  wr  for  any  values  of  k  vice  setting  wr= [44700 
rad/sec  -->  (7120  Hz)  first  natural  freq] 
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wt= ( (kl+2*k3) /m) A . 5;  %calculates  wr  for  any  values  of  k  vice  setting 
wt=[195000  rad/sec  -->  (31000  Hz)  second  natural  freq] 

er= (cl ) / (wr*m) ; %calculates  er  for  any  values  of  k  vice  setting  er=[0.89 
-->  damping  ratio  of  the  first  resonant  mode] 

et= (cl+2*c3) / (wt*m) ; %calculates  et  for  any  values  of  k  vice  setting 
et=[1.24  -->  damping  ratio  of  the  second  resonant  mode] 

O, 

o 

o, 

o 

%  Eventually  want  to  allow  user  to  easily  change  these 
angle  =  0; 

%  angle=30 ; %degrees  -->  angle  of  the  incident  sound  wave 
theta= (angle/180) *pi;  %rad-->  angle  of  the  incident  sound  wave 
w=2*pi*10000;  %rad/sec  -->  10000  Hz)  Excitation  frequency 
P=l;%Pa  -->  Pressure  Amplitude  of  the  Sound  Wave 

o, 

o 

%Calculate  the  values  that  do  not  change  with  time 


tao=d*sin (theta) /344 ;  %sec  --> 

delay  time 

due  to 

angle 

denomr= 
Ar  calc 

(wr' 

'2-w' 

' 2 ) 

A2+ (2*wr*er*w) 

A2 ; %  inside 

sqrt 

of 

Denominator 

of 

the 

denomt= 

(wt' 

'2-w ' 

'2) 

A2+ (2*wt*et*w) 

A2 ; %  inside 

sqrt 

of 

Denominator 

of 

the 

At  calc 

Ar=P*s/m* (sin (w*tao/2) /sqrt (denomr) ) ; %Amplitude  of  the  response  of  the 
rocking  mode 

At=P*s/m* (cos (w*tao/2 ) /sqrt (denomt) ); %Amplitude  of  the  response  of  the 
translating  mode 

Phir=-atan ( (2*wr*er*w) / (wr A2-wA2 ) ) ; %  Phase  constant  rocking 
Phit=-atan ( (2*wt*et*w) / (wtA2-wA2) ) ; %  Phase  constant  translating 

O, 

o 

o, 

o 

%  Set  our  time  steps 

te=500e-6;  %time  end  in  seconds  [This  value  is  arbitrary  and  can  be 
changed] 

t=0 : 2 . 5e-6 : te; %time  from  0  to  time  end 

o, 

o 

%Calculate  the  displacements  Y1 (t)  and  Yr(t)  [Paper  uses  X,  but  to 
%translate  to  the  normal  X  and  Y  plane,  I  am  calling  it  Y 
for  count  =  1 : 1 : (te/2 . 5e-6) +1 

Y1 (count) =8e6* (At* sin (w. *t (count) +Phit)  + 

Ar*cos (w. *t (count) +Phir) ) ;  %left  membrane  displacement  as  a  function  of 
time  in  micrometers 

XI (count) =-(((. 5)A2-Y1 (count)  .  A2)  .  A  .  5)  ; 

Yr (count) =8e6* (At*sin (w. *t (count) +Phit)  - 
Ar*cos (w. *t (count) +Phir) ) ;  %right  membrane  displacement  as  a  function 
of  time  in  micrometers 

Xr (count) = ( ( .5) A2-Yr (count) . A2) . A . 5; 

end; 

xarrow= . 3* sin (theta) ; 
yarrow= . 3*cos (theta)  +  .2; 

%  set  up  for  avi  movie: 

scrsz  =  get ( 0 , ' ScreenSize ' ) ;  %Gets  dot  length  and  width  of  screen 

screensize=2  %Size  can  be  adjusted  by  changing  this  variable.  The 

larger  the  #  the  smaller  the  screen.  1  =  entire  screen. 

fig  =  f igure (' Position ',[ 1  1  scrsz (3) /screensize  scrsz (4) /screensize] ) ; 

%sets  the  size  of  the  figure  to  fill  up  most  of  the  screen. 

set ( fig, ' DoubleBuf f er ' , ' on ' ) ; 
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mov  =  avif ile ( ' f lymovie . avi ' ,  'COMPRESSION',  'none') 

%  draw  the  first  position  of  the  wings,  represented  as  3  points  across 
%hwing=plot ( [XI (1)  0  Xr ( 1 ) ] , [ Y1 ( 1 )  0 
Yr ( 1 )  ]  ,  ' LineWidth ' , 4 ,  ' erasemode ' ,  ' xor ' ) 
hleft=plot ( [XI (1)  0 ]  ,  [ Y1 ( 1 )  0 

]  ,  ' LineWidth '  ,  8 ,  ' erasemode ',  'xor',  'color',  '  b  '  ) 
hold  on 

hright=plot ( [  0  Xr(l)],[  0 

Yr ( 1 ) ]  ,  ' LineWidth ' , 8 ,  ' erasemode ',  'xor',  'color',  ' r ' ) 
hold  on 

harrow=plot ( [  0  xarrow] , [  0.2 

yarrow] , ' LineWidth ' , 4 , ' erasemode ' , ' xor ' , ' color ' , ' k ' ) 

title  ('Plot  of  Membrane  Displacement  for  10000kHz,  1 

Pa ' , ' Font Size ' , 14 ) ; 

axis ([-.6  .6  -.5  .5]); 

xlabel ( ' Length  (mm) ' ) ; 

ylabel ( 'Height/100  (nm)  ' )  ; 

text(-.55,  0 . 4 ,[' Direction  of  sound  ' ] , ' FontSize ' , 16) 
text(0.4,  0.4, ['  Deg  '],' FontSize ', 16) 

angletext  =  text ( 0 . 35 , 0 . 4 , num2str (angle) ,' FontSize ', 1 6) 

%  Loop  through  wing  movement  for  each  wing: 
for  i=2 : length (XI ) 

%change  data  given  to  handle  hwing: 

%  set (hwing,  ' xdata ' ,  [XI ( i )  0  Xr ( i ) ] ,  ' ydata ' ,  [Y1  (i)  0  Yr  (i) ] ) ; 
set(hleft,  ' xdata ' ,  [XI (i)  0] ,  ' ydata ' ,  [Y1 (i)  0] ) ; 
set (hright, ' xdata ', [0  Xr ( i )],' ydata ',[ 0  Yr(i)]); 
set (harrow,  ' xdata ',[ 0  xarrow] ,' ydata ',[. 2  yarrow]); 

F  =  getframe (gca) ; 
mov  =  addf rame (mov, F) ; 
drawnow 

%  Slow  it  down  with  a  loop: 

%  for  slower  =  1:50000 

%  muchslower  =  slower  *  exp(2); 

%  end 
end 


for  angle  =  30:30:60 

%  angle=30 ; %degrees  -->  angle  of  the  incident  sound  wave 
theta= (angle/180) *pi;  %rad-->  angle  of  the  incident  sound  wave 
w=2*pi*10000;  %rad/sec  -->  10000  Hz)  Excitation  frequency 
P=l;%Pa  -->  Pressure  Amplitude  of  the  Sound  Wave 

o, 

o 

%Calculate  the  values  that  do  not  change  with  time 

tao=d*sin (theta) /344 ;  %sec  -->  delay  time  due  to  angle 

denomr= (wr A2-wA2 ) A2+ (2*wr*er*w) A2 ; %  inside  sqrt  of  Denominator  of  the 

Ar  calc 

denomt= (wt A2-wA2 ) A2+ (2*wt*et*w) A2 ; %  inside  sqrt  of  Denominator  of  the 
At  calc 
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Ar=P*s/m* (sin (w*tao/2) /sqrt (denomr) ) ; %Amplitude  of  the  response  of  the 
rocking  mode 

At=P*s/m* (cos (w*tao/2 ) /sqrt (denomt) ); %Amplitude  of  the  response  of  the 
translating  mode 

Phir=-atan ( (2*wr*er*w) / (wrA2-wA2) ) ; %  Phase  constant  rocking 
Phit=-atan ( (2*wt*et*w) / (wt A2-wA2 ) ) ; %  Phase  constant  translating 

O, 

o 

xarrow= . 3* sin (theta) ; 
yarrow= . 3*cos (theta)  +  .2; 

O, 

o 

%  Set  our  time  steps 

te=500e-6;  %time  end  in  seconds  [This  value  is  arbitrary  and  can  be 
changed] 

t=0 : 2 . 5e-6 : te; %time  from  0  to  time  end 

o, 

o 

%Calculate  the  displacements  Y1 (t)  and  Yr(t)  [Paper  uses  X,  but  to 
%translate  to  the  normal  X  and  Y  plane,  I  am  calling  it  Y 
for  count  =  1 : 1 : (te/2 . 5e-6) +1 

Y1 (count) =8e6* (At* sin (w. *t (count) +Phit)  + 

Ar*cos (w. *t (count) +Phir) ) ;  %left  membrane  displacement  as  a  function  of 
time  in  micrometers 

XI (count) =- ( ( ( .5) A2-Y1 (count)  .  A2)  . A . 5)  ; 

Yr (count) =8e6* (At*sin (w. *t (count) +Phit)  - 
Ar*cos (w. *t (count) +Phir) ) ;  %right  membrane  displacement  as  a  function 
of  time  in  micrometers 

Xr (count) = ( ( . 5) A2-Yr (count) . A2) . A . 5; 

end; 


set (angle text, 'string ' , num2str (angle) ) 

%  Loop  through  wing  movement  for  each  wing: 
for  i=2 : length (XI) 

%change  data  given  to  handle  hwing: 

%  set (hwing,  ' xdata '  ,  [XI ( i )  0  Xr ( i ) ] ,  ' ydata ' ,  [Y1 (i)  0  Yr  (i) ] )  ; 
set (hleft,  ' xdata ' ,  [XI (i)  0] ,  ' ydata ' ,  [Y1 (i)  0] )  ; 
set (hright,  ' xdata ' ,  [0  Xr (i) ] ,  ' ydata ' ,  [0  Yr (i) ] ) ; 
set (harrow,  ' xdata ',[ 0  xarrow]  ,  ' ydata '  ,  [ 0 . 2  yarrow]); 

F  =  getframe (gca) ; 
mov  =  addf rame (mov, F) ; 

drawnow 


%  Slow  it  down  with  a  loop: 


q, 

o 

for  slower  =  1:50000 

q, 

o 

muchslower  =  slower 

q, 

o 

end 

end 

end 

mov  = 

close (mov) ; 
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APPENDIX  B  -  DR.  KIM’S  ANSYS  CODE 


Finish 

/clear, nostart 
/prep7 

*DIM,A,array,100,3 

*CREATE,ansuitmp 

A/READ, A(1,1),data1,m,,JIK, 3, 100 

(3F10.4) 

*END 

/INPUT, ansuitmp 


et, 3, solid  45 
ET,  1,136,1  ! 
ET,  2,138,1 
ET,  4,136,1  !0 
et,5,200,6 


lelement  for  membrane 
!  4-node  option,  High  Knudsen  Number 

!  Circular  hole  option,  High  Knudsen  Number 
!  4-node  option,  High  Knudsen  Number 


igeometric  parameters 

L=305  !690  !  !390  !um,  length  of  mem 

W=305  !um,  half  width  of  mem 

L1=75  137.5  !75  !um,  length  of  torsion  beam 

W1=7  !um,  half  width  of  torsion  beam 

t=3.5  !um,  thickness  of  Si 

tau=0.5+0.75 


parsav,1  ,para,txt  !  Save  all  parameter  data 

off  1=5 


!sJ1  =  1050 
is  12=600 


!  Plate  hole  location  to  the  length  direction 

!  Plate  hole  location  to  the  width  direction 


c_r=43  !60  !  Hole  radius 

d_el=2  !+1  !  Gap+alpha  to  adjust  dmp/3 

d_el2=2  !+1  !gap2+alpha  to  adjust  dmp/3 

d_l=30 

d_v=30 

L2=L-(c_r+d_l-off1 )  !um,  length  of  mem  for  hole 


hnumj=3  !4  ihalf  hole  number  to  lateral, x 

hnum_v=3  !4  Ihalf  hole  number  to  vertical, y 

s_l1=(2*c_r+d_l)*hnumj  !  %  Plate  hole  location  to  the  length  direction 
sJ2=(2*c_r+d_v)*hnum_v  ! 


pamb=0.1 

visc=18.3e-12 

pref=0.1 

mfp=64e-3 

Knud=mfp/d_el 


!  ambient  pressure  (MPa) 
!  viscosity  kg/(um)(s) 

!  Reference  pressure  (MPa) 
!  mean  free  path  (urn) 

!  Knudsen  number 
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!!!!!!!!!!!!Grating 

L_g=300  !200  llength  of  grating 
W_g=7  Iwidth  of  grating 
W_d=7  ! 

W_h=W_g/2+W_d  lhalf  width  of  hole 
n_g=8  lhalf  number  of  grating 


mp,visc,1  ,visc 
mp, vise, 2, vise 
mp, vise, 4, vise 
mp,ex,3,130e3 
mp,nuxy,3,0.28 
mp,dens,3,2.32e-15 
mp,damp,5,7.182955*e-5 


!  Dynamic  viscosity  gap 

!  Dynamic  viscosity  holes 
!  Dynamic  viscosity  gap 
Imodulus  of  Si 
Ipoissons  ratio  Si 
Idensity  of  Si 

!  Material  damping  (from  squeeze  film  results) 


r,1,d_el,,,pamb 

rmore,pref,mfp 


!  Real  constants  -  gap 


r,2,c_r,,,pamb 

rmore,pref,mfp 


!  Real  constants  -  hole 


r,4,d_el2,,,pamb 

rmore,pref,mfp 


!  Real  constants  -  gap 


mmmmiBuiLD  modeumm! 

Icreate  membrane 
!biock,0,L,-W,W,0,t 
!biock,0,-L,-W,W,0,t 
!block,-W1  ,W1  ,W,W+L1 ,0,t 
!block,-W1  ,W1  ,-W,-W-L1 ,0,t 


!vsel,all 

!vglue,all 

Ivplot 


rectng,-L,L,-W,W 
rectng,c_r+d_l,3*c_r+d_l,c_r,-c_r 
agen,hnum_l-1 ,2,,,sJ1/hnum_l 
rectng,-c_r-d_l,-3*c_r-dj,c_r,-c_r 
agen,hnumJ-1,hnum_l+1,,,-sJ1/hnum_l 
agen,hnum_v,2,1+2*(hnum_l-1),,,s_l2/hnum_v 
agen,hnum_v,2,1+2*(hnum_l-1),,,-s_l2/hnum_v 
ASBA,  1,  all 


K,5001  ,L-L2,w,0 

K,5002,L-L2,-w,0 

K,5003,-L+L2,w,0 

K, 5004,-L+L2,-w,0 

L,  5001, 5002, 

L, 5003, 5004, 

lnum=4*(2*hnumJ-2)*(2*hnum_v-1)+4  Mine  number 
anum=1+(2*hnum_l-2)*(2*hnum_v-1)  larea  number 
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LSEL,s,LINE„lnum+1,lnum+2 
ASBL,anum+1 ,  all,  .delete,  delete 
!ASBL, 86, all,, delete, delete 

rectng,-W1  ,W1  ,W,W+L1 
rectng,-W1  ,W1  ,-W,-W-L1 


Igrating 

*do,jj,1,n_g 

rectng,L,L+L_g,W_h*(1+2*(jj-1))+W_g*(jj-1),W_h*(1+2*(jj-1))+W_g*jj 
rectng,L,L+L_g,-W_h*(1+2*(jj-1))-W_g*(jj-1),-W_h*(1+2*(jj-1))-W_g*jj 
rectng,-L,-L-L_g,W_h*(1  +2*(jj-1  ))+W_g*(jj-1  ),W_h*(1  +2*(jj-1  ))+W_g*jj 
rectng,-L,-L-L_g,-W_h*(1+2*(jj-1))-W_g*(jj-1),-W_h*(1+2*(jj-1))-W_g*jj 
*enddo 

aglue,all 

TYPE,  4 
MAT,  4 
esize,10 
!smrtsize,1 
AMESH,  39,40 


TYPE,  1 
MAT,  1 
esize,20  110 
!smrtsize,1 
IAMESH,  1,2,1 
AMESH,  4,5,1 
AMESH,  38 


TYPE,  5 
MAT,  5 
smrtsize,10  !5 
AMESH,  6,37 

!  Begin  Hole  generation 
*do,j,1,(2*hnum_v-1) 

*do,i,1,(hnum_l-1) 
k=i+(j-1)*(hnum_l-1) 
nsel  ,all 

*GET,  numb,  node,  ,  num,  max  !  Create  nodes  for  link  elements 

N,  numb+1,(i+0)*s_l1/hnum_l,-s_l2+j*s_l2/hnum_v, 

N,  numb+2,(i+0)*sJ1/hnum_l,-s_l2+j*s_l2/hnum_v,  t 

TYPE, 2 
MAT,  2 
REAL, 2 
NSEL,  all 

E,  numb+1,  numb+2  !  Define  2-D  link  element 
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ESEL,  s,  type, ,  1 

ESEL,a,type,,4 

NSLE,s,1 

local, 1 1 ,0,(i+0)*s_l1/hnumJ,-s_l2+j*s_l2/hnum_v 
csys,1 1 

NSEL.r,  loc,  x,  -c_r,c_r  !  Select  all  nodes  on  the  hole  circumference 
NSEL.r,  loc,  y,  -c_r,c_r  !  Select  all  nodes  on  the  hole  circumference 
NSEL,a,  node,  ,numb+1 

*GET,  next,  node,  ,  num,  min 
CP,  k,  pres,  numb+1,  next 
nsel,u,node,  , numb+1 
nsel,u,node,  .next 

CP,  k,  pres, all  ICoupled  DOF  set  for  constant  pressure 

csys.O 

*enddo 

*enddo 


*do,j,1,(2*hnum_v-1) 

*do,i,1,(hnum_l-1) 
kk=k+i+(j-1  )*(hnum_l-1 ) 
nsel.all 

*GET,  numb,  node,  ,  num,  max  !  Create  nodes  for  link  elements 

N,  numb+1, -(i+0)*sJ1/hnumJ,-s_l2+j*sJ2/hnum_v, 

N,  numb+2,-(i+0)*s_l1/hnumJ,-s_l2+j*s_l2/hnum_v,  t 

TYPE, 2 

MAT,  2 

REAL, 2 

NSEL,  all 

E,  numb+1,  numb+2  !  Define  2-D  link  element 
ESEL,  s,  type, ,  1 
ESEL, a, type, ,4 
NSLE,s,1 

local, 1 1 ,0,-(i+0)*sJ1/hnum_l,-s_l2+j*s_l2/hnum_v 
csys,1 1 

NSEL.r,  loc,  x,  -c_r,c_r  !  Select  all  nodes  on  the  hole  circumference 
NSEL.r,  loc,  y,  -c_r,c_r  !  Select  all  nodes  on  the  hole  circumference 
NSEL, a,  node,  , numb+1 

*GET,  next,  node,  ,  num,  min 
CP,  kk,  pres,  numb+1,  next 
nsel,u,node,  , numb+1 
nsel,u,node,  .next 

CP,  kk,  pres, all  ICoupled  DOF  set  for  constant  pressure 

csys.O 

*enddo 

*enddo 

!  End  hole  generation 
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!  Extrude  structural  domain 


esize,,2 
type,  3 
mat, 3 
real, 3 
vext,all . t 

esize,,2  !1 
type,  3 
mat, 3 
real, 3 

vext,208 . tau  !  Extrude  structural  domain 

!vext,1,„,,tau  !  Extrude  structural  domain 

!vext,43,,,,,tau  !  Extrude  structural  domain 

nsel,s,loc,x,L 

nsel,a,loc,x,-L 

nsel,a,loc,y,-W-L1 

nsel,a,loc,y,W+L1 

nsel,r,loc,z,-1e-9,1e-9 

d, all, pres,  !  Fix  pressure  at  outer  plate  boundary 

nsel,all 

nsel,s,loc,y,W 
nsel,a,loc,y,-W 
nsel,u,loc,x,-W1  ,W1 
nsel,r,loc,z,-1e-9,1e-9 

d, all, pres,  !  Fix  pressure  at  outer  plate  boundary 

nsel,all 

nsel,s,loc,x,W1 

nsel,a,loc,x,-W1 

nsel,u,loc,y,-W,W 

nsel,r,loc,z,-1e-9,1e-9 

d, all, pres,  !  Fix  pressure  at  outer  plate  boundary 

nsel,all 

esel,s,type,,2 

nsle,s,1 

nsel,r,loc,z,t 

d, all, pres, 0  !  P=0  at  top  of  plate 

!dlist,all 

esel,s,type,,1 

esel,a,type,,4 

nsle,s,1 

!nsle,a,4 

!nsle,s,all 

nsel,u,cp,,1 ,2*(hnumJ-1  )*(2*hnum_v-1 ) 

cm,FLUN,node 

allsel 

Imembrane  boundary  condtions 

asel,s,loc,y,W+L1-0.01  ,W+L1+0.01  Iselect  area  parallel  to  x  axis 

asel,a,loc,y,-W-L1-0.01,-W-L1+0.01  Iselect  area  parallel  to  x  axis 

nsla,s,1 
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d,ail,ux,0 
d,all,uy,0 
d,all,uz,0 
allsel 

finish 

/solu 

!antyp,harm  !  Full  Harmonic  analysis 

!harfrq,4000 

antype, modal 
modopt,lanb,2 
eqslv, sparse 
mxpand,2 
solve 
fini 

/postl 

RMFLVEC  !  Extract  eigenvectors 

fini 

/solu 
frl  =  1 

fr2=60000 

DMPEXT,  1 , 1  ,d_damp,fr1  ,fr2, 1 0 
DMPEXT,2,2,d_damp,fr1  ,fr2, 1 0 
!DMPEXT,3,3,d_damp,fr1,fr2,10 

/POST1 
SET, LIST, 999 

/solu 

abextract,1,2  !  Extract  damping  ratios  and  Rayleigh  constants 
finish 


!  Modal  analysis 

!  Extract  lowest  two  eigenmodes 
!  Expand  lowest  two  eigenmodes 
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APPENDIX  C  -  PROOF  OF  FWV2M 


nix  +  bx  +  kx  =  Fo  sin  cot 


...  2  F  . 

x  +  /x  +  a>~x  =  —sin  cot 
m 

y  =  b  /  m  coo2  =  k!  m 


F  /  m 


[( coo2-co2)2+(rco )2]2 

_ FJjmo)2) _ 

1 

[(i-V  2)2  +  (r(0/  2)2Y 

/CO  /  CO 


a=ay  r=  y 

/  CO  /  CO 


=>  XQ 


F0!(mco o2) 


[(l-Q2)2+(r2Q2)2]2 


V  -cox 


r  co 

— — (=Q) 

{mcoo)  coo 


V0=' 


[(i-n2)2+(r2Q2)2]2 

F0/(mco0 ) 


[(1ir)2+(r’)2]5 

v  (max)  =  EJXaeO  =>Q  =  l,a»  =  ft> 


r 


v0 (max)  2  |  (1-Q2.2  ' 

Q 

v„  1 


v0  (max) 

rn 
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So  if  — — —  =  ^then: 

vo  (max)  V  2 

va  _  1  _  1 
v0(max)  n  |  (1~QN21i  ^ 

1  1  TQ  1 

So  denominator  must  equal  2,  thus v  pq  ’  .  Therefore,  take  square 


root  of  both  sides: 


,(l-Q2)  =  ±rQ 


Q2+rQ-l  =  0 


Negative  gamma  is  not  possible 


Q2  -Qj  =  r 

^.  =  J-^FW42M  =  y  =  2B 
co„  C0o 
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